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Abstract 

The  Self-Induction  Theory  of  Vortex  Breakdown 
by  Charles  B.  Cain 

Chairperson  of  the  Supervisory  Committee 

Professor  Mitsuru  Kurosaka 
Department  of  Aeronautics  and  Astronautics 

A  new  study  of  the  transient  stages  leading  to  the  formation  of  vortex  breakdown 
shows  that  vortex  breakdown  is  initiated  by  a  negative  vorticity  gradient  that  causes  an 
inviscid  self-induction  feedback  mechamsm  resulting  in  steady  state  vortex  breakdown. 
We  call  this  the  self-induction  theory  of  vortex  breakdown.  The  vortex  filament  method 
captures  the  evolution  of  this  transient  formation  of  vortex  breakdown.  An  axial  vorticity 
gradient  is  introduced  into  the  vortex  tube  by  changing  the  circulation  along  the  tube. 
Thereafter,  the  self-induction  process  starts  on  its  own  as  the  axial  vorticity  induces 
azimuthal  velocity,  which  in  turn  tilts  the  vorticity  vector  in  the  azimuthal  direction.  Due 
to  the  gradient  in  azimuthal  vorticity  caused  by  the  increase  in  circulation,  the  vortex  tube 
radially  expands  and  the  vortex  filaments  contract  in  an  action  we  call  pile-up.  This  is 
followed  by  a  sign  switch  in  the  azimuthal  vorticity  caused  by  the  region  downstream  of 
the  vorticity  gradient  rotating  slower  than  the  upstream  region.  These  actions  proceed 
together  until  they  form  what  we  call  the  turning  point  where  the  vortex  filaments  turn 
inward  on  themselves  causing  a  sign  switch  in  the  axial  vorticity.  Vorticity  and  velocity 
data  produced  from  this  simulation  compare  well  to  experimental  data.  In  conjunction 
with  the  computer  simulation,  we  have  verified  these  results  experimentally  with  a  delta 
wing  model  in  a  water  tunnel  using  dye  flow  visualization,  laser-induced  fluorescence, 
and  particle  image  velocimetry.  These  results,  combined  with  comparisons  with  previous 
experiments  agree  with  one  another  and  support  the  self-induction  theory  of  vortex 
breakdown. 
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1.  Introduction  to  Vortex  Flows 

Vortex  flows  are  regions  of  rotating  flow  in  air,  water,  or  any  fluid.  Vortex  flows 
are  common  in  nature.  We  see  them  in  weather  phenomena  like  tornadoes  and 
hurricanes.  We  might  also  notice  a  vortex  that  develops  when  we  drain  a  sink  or  a  tub. 
Figure  1-1  shows  a  vortex  flow  in  an  image  of  Hurricane  Fran  taken  by  NASA-GSFC, 
using  the  NOAA  GOES-8  imaging  satellite. 


Figure  1-1,  Vortex  flows  in  nature  {ref.  1) 

Vortex  flows  are  also  commonly  found  around  aircraft,  especially  aircraft  with 
significant  leading  edge  sweep  or  leading  edge  strakes.  Figure  1-2  displays  two  images 
of  an  F/A-1 8  fighter  aircraft  with  a  vortex  generated  at  the  root  of  its  leading  edge  strake. 


Figure  1-2,  Vortex  flow  over  F/A-18  (ref.  2) 

At  these  high  angles  of  attack,  the  leading  edge  strokes  generate  two  strong 
vortices  on  the  suction  surface  of  the  wing.  In  the  first  pictures,  the  vortices  are  made 
visible  by  injecting  smoke  into  the  flow.  In  the  second  picture,  water  vapor  condenses  to 
make  the  vortices  visible. 

A  vortex  flow,  like  that  on  the  suction  surface  of  an  F/A-18  wing  generates  a 
strong  component  of  axial  velocity  along  the  centerline  of  the  vortex.  This  is  beneficial 
to  the  aircraft’s  performance  at  high  angles  of  attack  because  the  higher  velocities 
decrease  pressure,  thereby  increasing  lift.  The  kinetic  energy  of  the  vortex  also  energizes 
the  flow  in  the  boundary  layer  to  help  delay  separation.  The  overall  result  of  the  leading 
edge  strake  and  its  associated  vortex  is  an  increase  in  the  aircraft’s  maximum  angle  of 
attack  and  its  maximum  lift  coefficient  (ref.  3).  Certainly,  there  is  a  penalty  to  be  paid  for 
this  increased  lift  in  terms  of  increased  drag,  but  most  modem  fighter  aircraft  have 
sufficient  thmst  to  overcome  this. 

At  large  angles  of  attack,  this  vortex  can  undergo  a  rapid  change  in  flow 
properties  called  vortex  teeakdown.  Figure  1-3  shows  an  example  of  vortex  breakdown 
over  a  Plexiglas  delta  wing  in  a  water  tunnel  (ref.4). 


Figure  1-3,  Vortex  breakdown  over  a  delta  wing 

Vortex  breakdown  is  usually  characterized  by  a  rapid  radial  expansion  of  the 
vortex,  downstream  turbulence,  and  backflow.  This  rapid  change  in  flow  properties  can 
affect  the  stability  of  the  aircraft  and  the  downstream  turbulence  can  cause  structural 
damage  through  high-cycle  fatigue.  This  structural  damage  is  especially  prevalent  on 
twin-tailed  aircraft  like  the  F/A-18,  where  the  vertical  tail  is  immediately  downstream  of 
the  vortex.  This  is  clear  in  the  pictures  of  the  F/A-18  in  Figure  1-2. 

The  purpose  of  this  research  is  to  understand  how  we  can  preserve  the  positive 
effects  of  the  vortex,  while  suppressing  the  negative  effects  of  vortex  breakdown.  We  do 
this  by  trying  to  understand  the  cause  of  vortex  breakdown.  Many  others  have 
approached  this  same  problem  and  yet,  as  stated  recently,  “a  comprehensive  theory  of 
vortex  breakdown  still  does  not  exist”  (ref.  5).  Similar  comments  were  made  in  1977 
when  another  person  said,  “the  embarrassing  number  of  different  theoretical  notions  has 
not,  it  must  be  admitted,  led  to  satisfactory  understanding  of  the  flows  observed”  (ref.  6). 
So  there  is  work  yet  to  be  done. 

It  is  often  easier  to  study  vortices  and  vortex  breakdown  in  controlled  laboratory 
settings.  In  a  laboratory,  vortices  are  often  generated  in  pipes,  or  in  water  and  wind 
tunnels  with  a  simple  delta  wing  model. 
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1.1  Vortex  flow  in  a  pipe 

You  need  two  things  to  generate  a  vortex  in  a  pipe — something  to  introduce 
rotation  into  the  flow,  and  something  to  send  that  rotating  flow  downstream.  Usually, 
experimenters  use  swirl  vanes,  a  swirl  manifold,  or  angled  jets  do  this.  Swirl  vanes  are 
nothing  more  than  small  wing  sections,  arranged  around  the  circumference  of  the  pipe. 
Usually,  you  can  change  their  angle  of  attack  to  change  the  amovmt  of  swirl  in  the  flow. 
The  downwash  generated  by  the  wing  sections  deflects  the  flow  to  make  it  swirl.  A 
swirling  manifold  has  spiraling  ducts  that  introduce  rotation  into  the  flow.  Angled  jets 
send  fluid  into  the  pipe  at  an  angle  to  the  walls  of  the  pipe  to  cause  the  flow  to  rotate. 

The  velocity  component  that  sends  the  flow  downstream  can  be  created  any 
number  of  ways.  One  way  is  to  use  gravity.  The  pipe  is  setup  vertically  with  a  tank  of 
water  about  it.  The  swirl  generators  are  place  at  the  junction  between  the  pipe  and  the 
tank,  with  some  type  of  flow  control  valve  in  between  the  tank  and  the  pipe.  The  valve  is 
opened,  and  the  flow  moves  into  the  pipe,  through  the  swirl  generators,  and  down  to  a 
test  section  where  the  vortex  can  be  observed. 

Figure  1-5  shows  the  geometry  of  a  swirl  vane  vortex  generator.  Other  ways  to 
generate  the  downstream  velocity  component  involve  using  pumps,  fans,  or  exploiting 
flow  mechanisms  like  Ekman  layer  pumping  where  a  rotating  bottom  in  a  closed  cylinder 
causes  flow  to  circulate. 


Figure  1-4,  Swirl  vane  vortex  generator  in  a  pipe  (ref.  7) 
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Figure  1-5,  Numerical  model  of  swirl  vane  vortex  generator  in  a  pipe  {ref.  7) 
Vortex  flows  in  pipes  are  very  useful  for  study  because  the  vortices  can  be  made 
very  symmetric.  This  allows  the  researcher  to  clearly  see  the  rotating  flow  without  the 
asymmetries  and  turbulence  that  can  be  present  in  flow  over  wings  in  a  wind  or  water 
tunnel.  This  is  especially  important  when  studying  vortex  breakdown,  where  one  would 
like  to  separate  the  turbulence  and  potential  asymmetries  of  the  breakdown  from  the 
turbulence  and  asymmetries  of  the  surrounding  flow.  Figure  1-6  presents  an  example  of 
a  vortex  and  its  breakdown  in  a  pipe.  The  apparatus  used  to  generate  this  vortex  is  the 
one  shown  in  Figure  1-4.  The  symmetry  possible  with  this  type  of  apparatus  is  clear. 


Figure  1-6,  Vortex  breakdown  in  a  pipe  {ref.  8) 
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1.2  Vortex  flow  over  a  delta  wing 

Vortex  flows  in  a  pipe  can  be  very  useful  to  the  researcher.  However,  the 
apparatus  used  is  very  specialized.  Also,  the  vortex  that  forms  over  real  aircraft  will  not 
be  as  symmetrical  as  the  one  in  the  pipe.  Furthermore,  it  can  be  difficult  to  collect  data 
from  pipe  flows,  due  to  the  confined  nature  of  the  pipe.  For  that  reason,  we  also  study 
vortex  flows  in  regular  wind  and  water  tunnels.  A  simple  flat  plate  delta  wing  like  those 
in  Figure  1-7,  can  he  used  to  generate  vortices. 


Figure  1-7,  Horizontally  moxmted  delta  wing  in  water  tunnel 

Flow  over  a  triangular  shaped  delta  wing  is  characterized  by  two  counter-rotating 
vortices.  The  vortices  form  as  the  flow  separates  from  the  sharp,  swept-back  leading 
edge  of  the  delta  wing  as  it  tries  to  move  from  the  high-pressure  lower  surface  of  the 
wing  to  the  low-pressure  upper  surfece  of  the  wing,  as  seen  in  Figure  1-8.  This 
separation  line  is  called  Si  in  the  figure.  Once  separated,  the  low  pressure  on  the  top 
surface  of  the  wing  reattaches  the  flow  at  line  Ai.  However,  the  flow  separates  again. 
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along  line  S2,  drawn  off  the  surface  of  the  wing  by  the  even  lower  pressures 
accompanying  the  fast  moving  flow  that  separated  from  Si. 


Figure  1-8,  Vortex  flow  over  a  delta  wing  (ref.  9) 

After  the  secondary  separation,  the  flow  turns  in  on  itself  and  wraps  up  into  the 
primary  vortex.  The  flow  turning  inward  on  itself  is  self-similar,  in  that  as  it  turns 
inward,  the  surrounding  flow  that  previously  turned  inward  will  sustain  the  flow  in  its 
inward  turning  direction.  It  continues  to  spiral  inward  at  a  smaller  and  smaller  radius, 
thus  forming  the  vortex.  This  second  separation  causes  a  secondary  rotating  flow  that 
reattaches  along  the  attachment  line  A2.  The  secondary  flow  then  separates  again, 
turning  back  inward  to  form  a  secondary  vortex.  While  the  flow  is  rotating,  there  is  an 
axial  (chordwise)  velocity  component.  That  component  is  undisturbed  by  the  vortices 
near  the  centerline  of  the  wing. 

As  we  will  describe  in  more  detail  in  the  next  section,  very  low  pressures 
characterize  the  center  of  these  vortices.  It  is  this  strong  low  pressure  that  produces  the 
aircraft  perfijrmance  benefits  of  a  vortex  flows.  However,  these  benefits  do  not  last 
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forever.  As  the  delta  wing  is  brought  to  larger  angles  of  attack  (angles  of  attack  greater 
than  10  degrees  or  so,  depending  on  leading  edge  sweep  among  other  factors),  these 
strong,  previously  stable  vortices  will  breakdown.  They  will  transit  from  a  state  of 
organized  rotation,  to  a  state  of  spiraling  turbulence,  characterized  by  rapid  radial 
expansion  and  backflow.  Look  again  at  Figure  1-6  on  page  5  to  see  how  dramatic  the 
transition  is.  If  we  can  prevent,  delay,  or  weaken  this  transition,  we  can  use  the  positive 
benefits  of  the  vortex  at  higher  angles  of  attack  for  increased  aircraft  performance. 

The  key  is  to  understand  why  the  vortex  breaks  down.  Then,  we  can  devise  ways 
to  prevent  that  breakdown.  We  will  introduce  new  ideas  as  to  why  vortices  breakdown 
occurs  in  later  sections.  First,  we  should  describe  the  vortex  flow  more  precisely.  With  a 
more  precise  understanding,  we  can  think  of  ways  to  model  the  flow.  If  we  can  produce 
a  model  of  the  vortex  that  over  time  seems  to  develop  the  properties  of  vortex 
breakdown,  we  can  look  at  the  model  for  clues  as  to  the  mechanism  for  that  breakdown 
and  see  if  it  applies  to  real  flows.  If  it  does,  we  will  know  where  to  intercede  to  stop  or 
weaken  the  vortex  breakdown  in  real  flows. 

1.3  Equations  of  motion 

As  with  many  flows  in  fluid  dynamics,  we  can  describe  the  motion  of  vortex  flows 
with  several  equations.  First,  let  us  look  at  vortex  flows  in  a  broader  sense.  Let  us  define 
vortex  flows  as  flows  with  large  velocity  components  rotating  in  a  circular  path.  This  fits 
the  examples  of  vortex  flows  previously  presented.  There  may  be  other  velocity 
components  in  the  flow,  but  that  velocity  around  centerline  characterizes  it.  That  velocity 
will  be  called  the  azimuthal  velocity  and  will  be  given  the  symbol  ue  to  signify  that  the 
direction  of  the  velocity  is  in  the  theta  direction  of  a  cylindrical  coordinate  system. 
Another  often  mentioned  term  is  vorticity.  Vorticity  (co)  is  a  measure  of  rotation  in  a 
flow.  It  is  defined  as  the  curl  of  the  velocity.  So  vorticity  is  vector  that  measures  rotating 
motion  in  a  flow.  It  has  units  of  (time)~^ .  Equations  (1-1)  through  (1-4)  show  the 
coinponents  of  vorticity  in  cylindrical  coordinates 
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Using  the  right  hand  rule,  you  can  see  the  sense  of  rotatioa  If  you  align  the 
thumb  of  your  right  hand  with  the  vorticity  vector  and  curl  your  fingers,  the  flow  will  be 
rotating  in  the  same  direction  you  are  curling  your  fingers.  Figure  1-9  presents  a  diagram 
of  a  rotating  fluid  element. 


A. 


Figure  1-9,  Diagram  of  vortex  flow 

Now,  lets  take  this  element  of  rotating  flow  and  extend  it  along  a  line,  as  in  Figure 
1-10.  The  local  vorticity  vector  will  be  aligned  with  this  line  and  the  flow  will  rotate 
tangentially  about  this  line  at  every  point.  We  can  call  this  a  vortex  line.  The  vorticity 
vector  is  always  aligned  with  a  vortex  line,  so  the  velocity  associated  with  that  vorticity  is 
always  known.  A  formula,  called  the  Biot-Savart  Law,  describes  the  velocity  vector 
associated  with  a  given  vorticity  vector.  This  law  will  be  described  in  detail  later,  as  it 
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allows  us  to  calculate  how  vortices  interact  with  each  other  and  how  they  interact  with 
themselves. 


Figure  1-10,  Vortex  line 

Next,  lets  collect  several  vortex  lines  together  in  a  group.  If  they  are  oriented  in 
the  same  direction,  we  can  call  this  group  of  vortex  lines  a  vortex  tube. 


dA 


Figure  1-11,  Vortex  Tube 


It  might  be  important  to  measure  the  strength  of  this  vortex  tube.  One  way  would 
be  to  add  up  vorticity  in  the  vortex  tube,  multiplying  each  vorticity  vector  by  the  cross- 
sectional  area  of  the  vortex  tube  that  it  occupies.  This  is  essentially  what  is  done.  We 
call  the  measure  of  vortex  strength  circulation.  It  is  given  the  symbol  T,  though  at  later 
points  we  may  refer  to  circulation  as  C.  Later,  we  will  distinguish  between  the 
circulation  of  the  entire  vortex  tube  and  the  circulation  in  an  individual  vortex  line. 
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This  measurement  of  vortex  strength,  F,  is  an  important  measure  of  vortex 
flows.  We  can  ftirther  describe  this  term  in  several  ways.  One  useful  way  to  think  of 
circulation  is  to  consider  it  in  terms  of  vorticity  flux.  It  can  be  defined  in  terms  of 
velocity  or  vorticity.  In  terms  of  velocity,  circulation  is  the  line  integral  of  velocity 
around  a  closed  contour,  C,  as  Equation  (1-5)  describes.  Equation  (1-6)  shows  how  this 
expression  can  be  put  in  terms  of  vorticity  using  Stokes  theorem.  This  expression  of 
circulation  in  terms  of  vorticity  shows  that  vorticity  can  be  thought  of  as  circulation  per 
unit  area  {ref.  10). 


II 
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Let  us  look  in  more  detail  at  the  way  a  real  vortex,  like  the  examples  presented 
earlier,  actually  behave.  Real  vortex  flows  typically  have  distinct  regions  of  rotation. 
Each  region  rotates  with  its  own  types  of  azimuthal  velocity  profiles.  Figure  1-12  shows 
the  velocity  and  vorticity  distribution  of  a  real  vortex.  A  real  vortex  is  often  idealized  by 
making  the  vorticity  profile  a  uniform  step  function  in  the  vortex  core.  Figure  1-13 
shows  this  idealization.  This  idealized  vortex  is  called  a  Rankine  vortex. 
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Figure  1-12  Real  vortex  velocity  and  vorticity  distribution  (ref  1.10) 


Figure  1-13  Rankine  vortex  velocity  and  vorticity  distribution  (ref  10) 


In  this  section  we  will  be  speaking  mostly  about  the  Rankine  vortex,  though  the 
concepts  also  apply  to  the  real  vortex.  The  only  change  is  the  vorticity  distribution.  We 
will  also  be  assuming  that  the  flow  is  ineompressible  and  inviseid.  Let  us  look  more 
closely  at  the  two  regions  of  rotation.  The  first  is  the  solid-body  rotation  region.  This 
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region  is  often  called  the  vortex  core.  This  region  starts  at  the  core  of  the  vortex  and 
extends  out  to  the  radius,  rc.  As  its  name  implies,  this  region  of  the  vortex  rotates  like  a 
solid  body,  with  an  azimuthal  velocity  that  is  proportional  to  the  radial  distance  to  the 
center  of  the  vortex.  Equation  (1-7)  shows  the  azimuthal  velocity  distribution,  where  ©o 
is  the  angular  velocity  of  the  solid  body  rotation  region.  If  the  only  velocity  component 
is  that  in  the  ayirmithal  direction,  we  can  see  that  the  vorticity  in  this  region  only  in  the  z- 
component,  as  Figure  1-9  suggests.  Circulation  in  this  region  is  defined  in  terms  of 
Equation  (1-5)  and  described  in  Equation  (1-9).  Because  the  circulation  is  non-zero,  this 
region  of  the  vortex  is  rotational.  In  rotational  flow,  the  vorticity,  ©  s  curl^)  s  V  x  m  ,  is 
non-zero. 
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The  second  region  of  rotation  is  the  irrotational  region.  As  it’s  name  implies,  this 
region  is  not  rotational.  Hence,  the  vorticity  is  zero  and  the  circulation  is  constant  in  this 
region.  Figure  1-12  and  Figure  1-13  show  the  rotational  and  irrotational  vorticity  profile. 
In  the  irrotational  region,  the  azimuthal  velocity  is  proportional  to  the  inverse  of  the 
radial  distance  from  the  center  of  the  vortex,  as  described  in  Equation  (1-10)  where  c  is 


some  constant. 
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If  the  only  velocity  component  is  the  azimuthal  velocity,  then  the  vorticity  is  zero 
and  the  circulation  is  constant,  as  Equations  (1-1 1)  and  (1-12)  show.  As  a  result  of 
Equation  (1-12)  we  can  now  find  the  proportionality  constant  in  Equation  (1-10).  As 
Equation  (1-13)  describes,  that  constant  is  the  circulation,  T.  These  results  indicate,  that, 
outside  of  the  vortex  core,  the  circulation  is  constant  and  independent  of  radial  distance 
from  the  vortex  center.  Another  outcome  of  this,  implied  by  the  definition  of  circulation 
in  Equation  (1-5),  is  that  in  close  loop  regions  that  do  not  contain  the  vortex  core,  the 
circulation  is  zero.  For  closed  loop  regions  that  do  contain  the  vortex  core  and  are 
outside  of  that  core,  the  circulation  is  equal  to  the  expression  in  Equation  (1-14),  where  rc 
is  the  radius  of  the  vortex  core,  which  is  a  constant  for  inviscid  vortex  flows. 
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The  outcome  that  circulation  is  constant  for  all  time  along  a  closed  path  for 
inviscid,  incompressible,  barotropic  flow  is  known  as  Kelvin’s  circulation  theorem. 


Mathematically,  this  is  e)q)ressed  as  =  0 .  Thus,  irrotational  flow  (the  flow  outside  of 

the  vortex  core)  will  remain  irrotational  for  all  time,  under  the  stated  assumptions. 

Another  important  theorem  to  vortex  flows  is  the  Helmholtz  vortex  theorem.  This 
theorem  applies  directly  to  the  vortex  lines  that  were  mentioned  previously.  Again, 
^ggiitning  inviscid,  incompressible,  barotropic  flow,  Helmhotz  shows  that  (1)  vortex  lines 
move  with  the  fluid  elements,  (2)  circulation  is  constant  along  vortex  line,  (3)  a  vortex 
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line  or  vortex  tube  must  end  at  a  boundary  or  form  a  loop  (it  cannot  end  in  a  fluid),  and 
(4)  the  strength  of  a  vortex  line  in  a  vortex  tube  remains  constant  for  all  time  (ref.  1 0). 

We  have  defined  how  rotating  flow  acts  and  how  we  can  define  that  rotating  flow 
by  vorticity  and  vortex  lines.  The  next  question  is,  how  does  this  vortex  line  interact  with 
its  environment?  How  does  in  interact  Avith  itself?  To  answer  this  question,  we  need  to 
look  at  the  governing  equations  that  show  how  vorticity  can  change  in  the  fluid.  This  is 
succinctly  defined  with  the  so-called  vorticity  equation.  Next,  we  want  to  know  exactly 
what  velocity  field  is  associated  with  a  given  vorticity  field.  The  Biot-Savart  law  defines 
this  velocity.  With  these  principles  in  mind,  we  can  look  at  models  for  the  vortex  flow. 

1.3.1  Vorticity  equation 

The  vorticity  equation  is  essentially  an  angular  form  of  the  momentum  equation. 

It  is  derived  from  that  momentum  equation  (1-15),  using  the  continuity  equation  (1-16), 
and  a  vorticity  solenoidality  condition  (1-17). 


div(u)  =  0 

(1-15) 

+  (m  .  v)/ =  —  V(p)  +  vV^M 

(1-16) 

Dtdt  P 

div(Gi)  =  div(V  X  m)  =  0 

(1-17) 

Assuming  incompressible  flow  and  no  body  force,  we  can  find  the  vorticity 
equation  by  taking  the  curl  of  the  momentum  (1-21)  equation  and  using  the  vector 
identities  in  Equations  (1-18),  (1-19),  and  (1-20)  (ref  11).  What  results  is  the  vorticity 
equation,  which  is  shown  in  Equation  (1-22). 
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curl{VA))^t) 

(1-18) 

div(curl(A))  =  V  •  (V  x  .3)  =  0 

(1-19) 

curl(A  X  5)  =  Adiv{B)  -  Bdiv{A)  +  {B  •  V)^  -  (A  •  V)B  =  0 

(1-20) 

curl 

+  (m  -  v)m  =  —  V(p)  +  vV"m 
[Dt  dt  ^  ^  p  J 

(1-21) 

^  +  fc .  vk  -  (®  •  v)m  +  vV'w 

Dt  dt  ^  ^  ^  ^ 

(1-22) 

There  are  several  interesting  things  in  Equation  (1-22).  First,  there  is  no  pressure 
term.  It  drops  out  because  the  curl  of  the  pressure  gradient  is  zero,  by  the  vector  identity 
(1-18).  If  the  fluid  is  barotropic,  then  the  center  of  pressure  always  acts  through  the 
center  of  mass  in  a  fluid  element.  Thus,  pressure  can  generate  no  torque  and  no  vorticity. 
If  we  had  included  conservative  body  forces  like  gravity,  they  would  have  dropped  out  of 
the  equation  for  the  same  reason. 

For  this  incompressible,  barotropic  flow,  vorticity  can  only  change  in  two  ways 
corresponding  to  the  two  terms  on  the  right  hand  side  of  (1-22).  Let  us  look  at  each  of 
these  terms.  First,  look  at  the  term  vV^ro .  The  v  term  is  the  kinematic  viscosity.  The 

solution  to - —  vV^co  is  like  the  solution  to  the  heat  equation.  Thus,  vorticity  diffuses 

Dt 

like  heat.  By  solving  this  differential  equation  for  a  two-dimensional  vortex  with  the 
vorticity  initially  concentrated  at  the  origin,  that  is  with  the  initial  condition 

(o(r,t)  =  ^  it  can  be  shown  that  the  core  radius  of  the  vortex  diffuses  according  to 

2nr 

=yl\.256-(4vt).  Inair,  v  =  1.45x10"^ and  in  water,v  =  1.12x10“^  — .  For 

s  ^ 

these  values,  it  would  take  that  imtially  concentrated  vortex  3.8  hrs  to  diffuse  to  a  core 
radius  of  1  meter  in  air.  In  water  it  would  take  49.4  hours  to  diffuse  to  that  radius.  If  the 
time  scale  of  the  motion  of  interest  are  considerably  shorter  than  time  scales  of  the  vortex 
diffusion,  it  is  acceptable  to  say  that  diffusion  is  small  and  its  effects  are  negligible. 
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The  second  term  in  this  vorticity  equation  is  (c5  •  v)/ .  This  is  called  the 
stretching  and  tilting  term,  and  it  describes  how  angular  momentum  is  conserved  when  a 
vortex  undergoes  these  processes.  Essentially,  if  a  vortex  tube  is  stretched  out,  its  radius 
will  decrease.  To  conserve  angular  momentum,  it  will  have  to  rotate  fester,  since  the 
moment  of  inertia  of  the  fluid  elements  decreases  as  it  is  stretched  out  to  conserve  mass. 
This  is  equivalent  to  the  familiar  example  of  ice  skaters  spinning  slowly  with  their  arms 
stretched  out  and  spinning  faster  when  they  draw  their  arms  to  their  bodies.  Likewise,  if 
a  vortex  is  tilted,  components  of  that  vorticity  will  translate  into  the  direction  the  vortex  is 
tilted.  This  means  that  the  vorticity  vector  is  at  all  times  proportional  to  the  fluid  line 
vector,  dr  .  Thus,  at  t  =  to ,  m  =  const.  •  ,  with  this  same  proportionahty  constant 

holding  for  any  t. 

Now  we  have  determined  that,  rmder  the  assumed  conditions,  once  introduced  to  a 
flow  field,  vorticity  only  can  change  by  tilting,  stretching,  and  difiusioa  Even  more,  we 
have  determined  that,  over  short  periods,  the  effects  of  diffiision  are  small.  Now  that  we 
have  vorticity  present,  and  know  how  it  can  change,  let  us  find  out  the  velocities  that  are 
associated  with  that  vorticity  field.  To  find  this,  we  turn  to  the  Biot-Savart  Law. 

1.3.2  Derivation  of  Biot-Savart  law 

The  Biot-Savart  law  gets  its  name  fi-om  a  magnetic  field  equation  that  describes 
the  magnetic  field  that  is  associated  with  electric  current  as  it  flows  through  a  'wire  {ref 
1 2).  In  one  of  those  wonderfiil  parallels  in  nature,  this  law  also  describes  the  velocity 
associated  with  a  given  vorticity  profile.  For  example,  if  we  specified  a  vortex  line  that 
formed  a  loop  with  a  given  radius  and  a  given  circulation,  we  could  find  that  velocity  that 
the  vortex  line  induces  at  the  center  of  that  loop,  or  at  any  point  around  that  loop.  In  a 
more  practieal  sense,  the  Biot-Savart  law  would  deseribe  how  the  wingtip  vortices  of  an 
aircraft  interact  'with  each  other,  as  in  Figure  1-14. 
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Dryden  Flight  Research  Center  ECN  7848  1977  Photo 

L-1011  Wi  ng  Vortex  Study.  NASA  photo 


Figure  1-14,  Wingtip  vortex  interactions  (ref.  13) 

The  Biot-Savart  law  is  very  important.  It  uniquely  describes  the  velocity  field 
associated  with  a  given  vorticity  field.  If  this  field  acts  linearly,  the  associated  velocity 
fields  of  several  different  vortex  lines  and  their  vorticity  fields  can  be  linearly  added  to 
each.  This  will  produce  an  overall  velocity  field  that,  while  in  sum  is  very  complicated, 
is  simply  the  linear  addition  of  many  smaller  parts  and  straightforward  calculations. 

Let  us  assume  that  vorticity,  a{x\y’,z'),  is  given  and  that  the  flow  is 
incompressible.  We  are  trying  to  find  the  relationship  between  the  given  ©at  point 
(x',y  ,z’)  and  the  velocity,  m  ,  at  a  different  point,  {x,y,z).  Since  the  flow  is 
incompressible,  we  know  (1-23)  is  true. 
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div(u)  =  0 

(1-23) 

u  =  curl  (A) 

(1-24) 

div{curl(A))  =  0 

(1-25) 

Then  vector  potential,  .i,  satisfying  (1-23)  exist.  Furthermore,  we  know  (1-25) 
by  the  vector  identity  (1-19).  It  should  be  noted  that  i  is  not  unique.  Ifweaddto^  any 
arbitrary  vector  Aq=  gf'£id(<ft) ,  (1-28)  will  result  by  the  vector  identity  (1-27)  (r^  JJ). 
Therefore  is  indeterminate.  As  a  result,  we  can  choose  Ji  as  we  like.  For  later 
convenience,  we  will  choose  such  that  (1-29)  is  true. 


A+Aq=  A  +  gradi^) 

(1-26) 

curl{grad{^y)  =  0 

(1-27) 

u  =  curl(A  +  grflt/((|>))  =  curl  (A)  +  curl{grad{^))  =  curl{A) 

(1-28) 

div(A)  =  0 

(1-29) 

Now,  we  can  substitute  (1-24)  into  the  definition  of  vorticity. 

as  in  (1-30).  Using 

the  vector  identity  (1-31)  (re/  ii),  we  find  (1-32).  Because  we  have  chosen  A  such  that 
it  satisfies  (1-29),  (1-32)  reduces  to  (1-33). 

®  =  curl(curl(Ajj 

(1-30) 

curl(curl(F  ))=  grad{div{F ))- 

(1-31) 

&  =  gradldiv{Ajj-V^A 

(1-32) 

G>  =  -V^A 

(1-33) 

If  we  split  (1-33)  into  its  vector  components,  we  get  the  expressions  in  (1-34), 
(1-35),  and  (1-36). 
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(1-34) 

1 

S 

11 

< 

> 

(1-35) 

(1-36) 

These  expressions  are  in  the  form  of  Poisson’s  equation  (1-37).  The  solution  of 
Poisson’s  equation  is  well  known  (1-38). 

V^i^{x,y,z)  =  f{x,y,z)  (1-37) 


The  integral  in  (1-38)  is  evaluated  over  the  entire  space.  Now,  let  us  define  two 
points  in  space,  as  in  Figure  1-15.  We  will  call  point  P  the  source  point  and  point  Q  the 
observation  point.  Thus,  we  are  going  to  observe  something  at  point  Q  whose  source  is  at 


point  P.  Let  us  also  define  the  vector  between  then  as  R . 


Figure  1-15  Source  and  observation  points 
If  we  apply  the  solution  to  Poisson’s  equation  to  the  expressions  in  (1-34),  (1-35), 
and  (1-36)  we  will  find  (1-39),  (1-40),  and  (1-41).  Note  that  the  denominator  of  (1-38)  is 

simply  the  magnitude  of  R  as  defined  in  Figure  1-15.  Also,  note  that  dV  =  dx  dy  dz  . 
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(1-39) 

(1-40) 

(1-41) 


Remembering  how  we  defined  A  in  (1-24),  we  can  break  that  up  into  its  vector 
components,  as  in  (1-42),  (1-43),  and  (1-44). 


u 


X 


Uy  =  curl{A)y  = 

u,  =CMr/(i),  = 


dA, 

dAy 

dz 

54 

dz 

dx 

dAy 

54 

dx 

dy 

(1-42) 


(1-43) 

(1-44) 


Now,  if  we  take  the  partial  derivatives  of  the  components  of  2 ,  as  defined  in 
(1-39)  through  (1-41),  and  substitute  those  partial  derivatives,  as  shown  below  in  (1-45) 
through  (1-50),  into  the  expressions  for  u  in  (1-42)  through  (1-44),  we  will  get  an 
expression  for  velocity  as  a  function  of  vorticity  and  distance. 
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-y^r'dV'  (1-45) 

(1-46) 

(1-47) 

^  Jf|0>;,(^'  >  /  .  ^')  -  ^’1^1 

dz  4n 

(1-48) 

Sx:  45T  •’•’•'  "  " 

(1-49) 

Collecting  all  these  ejqiressions  together,  and  recognizing  the  vector  components 
of  rox  ^  we  get  the  following  e>q>ressions,  (1-51)  to  (1-53),  for  m  as  a  function  of  ®  and 
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vortex.  If  we  are  interested  in  the  velocity  associated  at  point  Q  from  the  vorticity  at 
point  P,  as  in  Figure  1-12,  and  point  Q  is  far  enough  away  from  point  P,  this  does  not 
present  a  problem.  However,  if  we  need  to  find  the  velocity  close  to  the  vorticity  at  point 
P,  we  need  to  modify  (1-54). 

Another  factor  that  (1-54)  does  not  take  into  account  is  vortex  stretching.  If  we 
define  the  vortex  core  radius  as  a  =  ,  and  designate  a  constant  vortex  core  volume, 


based  on  an  initial  distribution  of  vorticity,  we  can  account  for  stretching  by  changing  a 
as  required  to  keep  the  core  volume  constant.  To  account  for  the  effects  of  the  vortex 
core  and  stretching,  we  modify  the  Biot-Savart  law  by  a  function,  s,  that  is  dependent  on 
R  and  the  vortex  core  diameter  (ref.  14).  Since  we  are  going  to  be  calculating  the 
velocities  induced  by  one  vortex  at  point  P  on  another  vortex  at  point  Q,  we  need  to 
account  for  the  diameter  of  both  vortices.  For  that  reason,  the  s-fiinction  is  dependent  on 
the  diameter  of  both  vortex  tubes.  The  results  of  these  modifications  are  in  (1-55)  and 
(1-56).  In  (1-56)  the  a-term  sets  the  distribution  of  vorticity  within  the  core.  Setting 
a  =  0.413  gives  a  Gaussian  distribution  of  vorticity  for  the  vortex  core  (ref  14). 


(1-55) 


|j4<T„a5)= 


(1-56) 


1.S.3  Physical  interpretation  of  equation  of  the  Biot-Savart  law 

With  the  expression  in  Equation  (1-55)  and  (1-56),  we  can  now  find  the  velocity 
at  the  observation  pomt  Q  that  is  associated  with  the  given  vorticity  at  the  source  point  P, 
taking  into  account  the  solid-body  rotation  vortex  core  and  the  effects  of  stretching.  This 
is  clear  mathematically,  but  what  does  in  mean  physically?  How  can  we  describe  motion 
in  various  type  of  fluids  without  any  viscous  or  pressure  terms? 
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To  answer  these  questions,  we  should  describe  just  what  the  Biot-Savart 
equation  is.  The  Biot-Savart  law  is  not  prescriptive  in  that  it  does  not  define  what  the 
velocity  field  will  be.  Rather,  it  is  descriptive.  It  describes  what  is.  It  simply  describes 
what  velocity  field  is  associated  with  a  given  vorticity  field  at  a  given  time.  We  could 
equally  describe  what  vorticity  is  associated  with  a  given  velocity  by  taking  the  curl  of 
the  velocity  field  (©  =  curliji)).  So  in  that  way,  the  Biot-Savart  law  is  a  translation  tool. 
We  can  describe  fluid  motion  in  the  language  of  velocity  or  in  the  language  of  vorticity. 
The  Biot-Savart  law  simply  tells  us  how  to  switch  between  the  two  languages. 

But  should  not  the  translation  depend  on  the  fluid  the  motion  is  in?  Do  we  not 
need  pressure  terms  to  move  fluid?  The  key  is  that  the  Biot-Savart  law  does  not  describe 
transient  conditions.  It  assumes  the  flow  is  steady.  When  steady,  the  viscosity  of  the 
fluid  mediums  and  the  pressures  have  done  their  job  by  bringing  the  fluid  into  a  steady 
state. 

However,  we  want  to  use  the  Biot-Savart  law  to  calculate  unsteady  time-evolution 
of  vortex  flows.  That  is,  we  want  to  take  a  given  vorticity  field  and  use  it  to  find  the 
associated  induced  velocities  at  points  around  it  (including  points  that  contain  vorticity). 
Then,  we  take  those  induced  velocities,  multiply  them  by  a  time  step,  and  calculate  how 
far  those  points  will  move.  We  then  move  them  by  that  amount  and  recalculate  the 
velocity  associated  with  the  new  distribution  of  vorticity.  We  should  he  able  to  use  the 
Biot-Savart  law  in  this  way  as  long  as  we  assume  that  the  associated  steady  induced 
velocity  field  is  formed  at  a  time  scale  smaller  than  the  time  step  used  to  move  the  flow 
solution  forward  in  time.  At  distances  on  the  order  of  the  diameter  of  most  vortex  tubes, 
this  should  be  true  in  most  fluid,  except  perhaps  fluids  with  extremely  low  viscosity,  like 
super-cooled  liquid  helium  or  the  like. 

Another  thing  we  should  assume  to  be  able  to  use  the  Biot-Savart  law  in  this  way 
is  infinite  (or  at  least  sufficiently  large)  signal  speed.  If  this  is  the  case,  there  time  lag  due 
to  signal  propagation  will  be  sufficiently  small  that  it  occurs  at  time  scales  much  smaller 
than  the  time  scale  used  to  advance  the  solution  in  time.  Again,  if  the  fluid  is 
incompressible  (or  nearly  so,  as  is  the  case  in  most  vortex  flows)  this  is  not  a  concern. 
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The  bottom  line  is  that  an  equation  is  only  as  good  as  its  assumption.  However, 
within  the  limits  of  those  assumptions,  an  equation  can  be  very  good  at  describing 
physical  effects,  even  in  very  simple  terms.  The  Biot-Savart  law  is  one  such  equation. 

As  long  as  we  meet  the  necessary  assumptions,  we  can  use  the  Biot-Savart  law  to 
describe  how  an  element  of  vorticity  interacts  with  the  fluid  surrounding  it,  to  include 
other  elements  of  vorticity.  We  will  use  this  property  to  break  a  vortex  tube  into  many 
vortex  elements  and  then  describe,  using  the  Biot-Savart  law,  how  those  vortex  elements 
interact  with  each  other.  This  is  known  as  the  vortex  filament  method.  If  we  can  use  this 
method,  with  an  appropriate  set  of  initial  conditions,  to  model  vortex  breakdown,  perhaps 
we  can  better  understand  the  mechanism  that  causes  that  vortex  breakdown.  From  this 
understanding  of  the  mechanism  that  causes  the  breakdown,  we  should  be  able  to  gain 
some  insight  as  to  how  we  can  interrupt  that  mechanism  to  delay,  prevent,  or  weaken 
vortex  breakdown  and  its  detrimental  effects. 
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2.  Vortex  filament  method 

The  vortex  filament  method  allows  us  to  look  at  vortex  flows  using  the  Biot- 
Savart  law,  as  was  previously  developed,  to  calculate  the  motion  of  those  flows.  If  the 
vortex  filament  method  captures  the  important  physics  of  the  flow,  it  will  yield  accurate 
results.  Conversely,  if  the  vortex  filament  method  yields  accurate  results,  we  can  know 
with  a  reasonable  degree  of  certainty,  that  the  important  physics  in  the  flow  are  limited  to 
those  present  in  the  vortex  filament  method. 

The  Biot-Savart  law  describes  the  induced  velocities  associated  with  a  given 
distribution  of  vorticity.  If  we  modify  that  law  to  include  the  induced  velocity  profile 
inside  a  vortex  core  and  the  efiects  of  vortex  stretching,  we  can  limit  our  discussion  of 
vortex  motions  that  are  calculated  using  the  vortex  filament  method  to  those  things 
associated  with  induced  velocities,  the  vortex  core  velocity  profile,  and  vortex  stretching. 
With  this  in  mind,  let  us  know  turn  to  how  we  can  use  the  vortex  filament  method  to 
model  a  vortex  tube. 

2.1  Vortex  tube  interacts  with  itself 

As  we  discussed  previously,  a  vortex  tube  is  composed  of  many  vortex  lines.  A  vortex 
line  is  the  same  thing  as  a  vortex  filament.  So  what  we  call  the  vortex  filament  method 
could  just  as  easily  be  called  the  vortex  line  method.  These  vortex  lines  are  aligned  with 
the  vorticity  vector  associated  with  the  fluid  motion,  as  in 

Figure  2-1. 
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Figure  2-1,  Vortex  Tube 

With  the  Biot-Savart  law,  we  should  be  able  to  calculate  the  velocity  induced  by 
the  vorticity  of  one  vortex  line  on  another  vortex  line.  With  this,  we  can  find  how  that 
other  vortex  line  will  move  as  a  result  of  this  induced  velocity.  By  the  Helmhotz  vortex 
theorem,  we  know  that  vorticity  will  move  with  the  motion  of  a  material  fluid  line.  This 
means  that  as  induced  velocities  cause  motion  to  fluid  elements,  they  in  turn  move  the 
vorticity  associated  with  those  fluid  elements.  For  this  reason,  we  can  calculate  the 
induced  velocities  associated  with  the  vortex  lines  of  a  vortex  tube,  multiply  those 
velocities  by  a  small  time  step,  and  then  calculate  how  fer  the  vortex  lines  will  move  in 
that  small  time  step.  We  do  this  again  and  again,  and  observe  the  vortex  tube  moving 
forward  in  time.  In  this  way,  the  vortex  lines  that  compose  a  vortex  tube  interact  with 
themselves  such  that  this  interaction  can  affect  the  shape  of  that  vortex  tube.  If  a  real 
vortex  tube  can  be  described  with  a  series  of  vortex  lines,  with  a  given  distribution  of 
vorticity,  then  the  vortex  filament  method  can  be  used  to  model  real  vortex  motion. 

2.2  Discretizing  the  vortex  tube 

To  use  the  vortex  filament  method,  we  need  to  discretized  a  vortex  tube  into  an 
adequate  number  of  vortex  filaments. 


Figure  2-2  Real  vortex  velocity  and  vorticity  distribution  (ref.  10) 
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A  vortex  tube  has  a  certain  cross-sectional  area,  defined  by  the  vortex  tube  core 
radius,  and  a  certain  distribution  of  vorticity  within  that  area.  Figure  2-2  shows  that 
distribution  for  a  real  vortex.  We  only  need  to  model  the  vortex  core  because  we  know 
that  almost  all  of  the  vorticity  is  contained  in  this  core.  Outside  of  the  core,  the  flow  is 
nearly  irrotational,  and  thus  contains  a  negligible  amount  of  vorticity.  A  vortex  line  with 
very  little  vorticity  will  not  induce  significant  velocities  and,  for  purposes  of 
discretization,  it  can  be  neglected. 

To  adequately  discretize  that  vortex  tube,  the  vortex  lines,  should,  in  total,  occupy 
the  same  cross  sectional  area  as  the  vortex  tube.  To  do  this,  each  vortex  line  is  given  its 
own  cross  sectional  area  and  its  own  vortex  core  radius.  Each  vortex  line  has  the  same 
Gaussian  distribution  of  vorticity  as  the  overall  vortex  tube.  It  is  the  vorticity  magnitude 
and  the  vortex  line  core  radius  of  each  vortex  line  that  is  set,  as  required,  so  that  in  sum, 
the  overall  vortex  tube  core  radius  and  vorticity  distribution  of  the  now  discretized  vortex 
tube  matches  that  of  the  original  vortex  tube.  These  vortex  lines  are  then  packed  densely 
enough  that  there  is  no  empty  space  inside  the  vortex  tube.  This  overlap  allows  the 
discretized  vortex  tube  to  maintain  its  original  smooth  distribution  of  vorticity  within  the 
original  vortex  tube.  Figure  2-3  shows  an  illustration  of  a  discretized  vortex  tube. 


Figure  2-3  Vortex  Tube  Discritization 

After  the  vortex  tube  is  discretized  into  many  overlapping  axial  vortex  lines,  it  is 
further  discretized  along  those  lines.  The  vortex  lines  are  broken  up  into  many  smaller 
vortex  filaments  segments.  This  allows  portions  of  the  vortex  line  to  be  moved  by 
induced  velocities  in  different  directions  along  its  axis.  The  character  of  the  original 
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vortex  lines  is  maintained  by  enforcing  the  condition  that  each  vortex  filament  segment 
is  connected  to  the  one  initially  upstream  and  downstream  of  itself.  This  preserves  the 
material  line.  We  will  refer  to  the  endpoints  of  each  vortex  filament  segment  as  a  nodes. 
The  induced  velocity  associated  with  all  other  vortex  filament  segments  is  calculated  at 
these  nodes.  This  total  induced  velocity  is  multiplied  by  the  time  step  so  that  the  location 
of  each  endpoint  node  of  the  vortex  filament  segment  can  move  in  time.  Figure  2-4 
shows  a  diagram  of  this. 
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Vortex  filament  segment 

Figure  2-4  Vortex  Filament  Segments 

This,  in  essence,  is  how  the  vortex  filament  method  works.  How  these  vortex 
lines  are  initially  setup  and  what  boundary  conditions  or  external  flow  conditions  are 
imposed  on  them  are  specific  to  each  individual  model.  Next  we  will  describe  the  how 
the  model  is  specifically  setup  for  this  research. 

2.3  Model  setup 

As  with  many  numerical  model,  our  vortex  filament  model  can  be  discussed  in 
terms  of  what  it  ideally  represents,  and  by  what  reasonable  compromises  we  make  to 
make  it  computationally  more  efficient.  The  explicit  calculation  of  induced  velocities  of 
one  vortex  filament  segment  on  every  other  vortex  filament  requires  (Xn  )  operations 
iref.  15),  so  we  hope  to  minimize  the  required  operations  whenever  possible.  We  shall 
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present  the  ideal  case  first,  and  then  discuss  what  changes  we  made  to  this  ideal 
representation  to  turn  the  ideal  representation  into  an  efficient  numerical  code. 


2.3.1  Ideal  case 

We  model  an  unconstrained  vortex  tube  that  is  being  fed  downstream  {z  >  O)  from 
a  vorticity  source,  located  at  z  =  0 ,  at  a  steady  rate  that  we  will  call  C/„  (z) .  The  vortex 
tube  convects  downstream  at  that  rate,  U„(z) .  All  vorticity  in  the  vortex  tube  is  initially 
aligned  in  the  axial  direction  (which  we  will  call  z).  As  it  is  initially  fed  out  of  the 
vorticity  source,  the  circulation  in  the  vortex  lines  is  set  such  that  there  is  overall  a 
Gaussian  distribution  of  this  axial  vorticity  in  the  vortex  tube.  Figure  2-5  shows  a 
diagram  of  this  vorticity  source.  The  circulation  in  the  vortex  tube  is  described  as  it  was 


previously,  C  =  jcurl(u)  dA  =  .  If  we  are  speaking  of  the  circulation  of  the  entire 

A  A 

vortex  tube,  we  will  call  that  circulation  C.  If  we  are  speaking  of  the  circulation  in  a 
particular  vortex  line,  we  will  call  the  circulation  T  . 


Figure  2-5  Axial  Vorticity  Source 

We  will  allow  the  circulation  of  the  vortex  tube,  as  it  comes  from  the  vorticity 
source  to  change  in  time.  Thus  at  one  time,  the  vorticity  source  may  produce  a  vortex 
tube  with  circulation  C  =  Cj ,  while  at  another  time  it  may  produce  a  vortex  tube  with  a 
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different  circulation  C  =  C2 .  At  this  point,  we  must  be  careful  not  to  violate  any  of  the 
Helmhotz  vortex  theorems. 

Let  us  restate  the  Helmholtz  vortex  theorem.  Assuming  inviscid,  incompressible, 
barotropic  flow,  Helmhotz  shows  that  (1)  vortex  lines  move  with  the  fluid  elements,  (2) 
circulation  is  constant  along  a  vortex  line,  (3)  a  vortex  line  or  vortex  tube  must  end  at  a 
boundary  or  form  a  loop,  it  cannot  end  in  a  fluid,  and  (4)  the  strength  of  a  vortex  line  in  a 
vortex  tube  remains  constant  for  all  time  {ref.  10). 

The  model  meets  the  first  part  by  moving  the  vortex  filament  segments  according 
to  the  velocity  induced  on  the  segment  endpoints  and  the  downstream  convection 
velocity,  U^^) .  The  second  part  is  met  until  we  change  the  circulation  that  is  produced 
by  the  vortex  source.  At  this  point,  since  we  cannot  change  the  vorticity  in  the  vortex 
lines  coming  from  the  vortex  source,  the  only  way  we  can  increase  the  circulation  is  to 
add  vortex  lines  to  the  vortex  tube.  We  can  place  these  vortex  lines  coincident  to  the 
existing  vortex  lines,  but  we  must  do  so  keeping  the  third  part  of  the  Helmholtz  vortex 
theorem  in  mind. 

Vortex  lines  and  vortex  tubes  cannot  end  in  a  fluid.  They  must  have  a  boundary,  as 
they  do  at  the  vorticity  source,  they  must  form  closed  loops,  or  they  must  exists 
downstream  to  z  =  00 .  What  our  model  will  do  is  insure  that  all  the  vortex  lines  that 
make  up  the  vortex  tube  form  closed  loops.  In  this  way,  the  vortex  tube  will  not  end  in  a 
fluid,  but  will,  itself  form  a  closed  loop  by  way  of  its  vortex  lines.  Let  us  make  this  clear 
with  the  illustration  in  Figure  2-6. 


Figure  2-6  Closed  Loop  Vortex  Lines 

In  this  example,  we  discretized  a  vortex  tube  into  four  vortex  lines.  Please  note 
that  at  this  time,  we  are  not  considering  the  time  evolution  of  the  vortex  tube.  If  it  were 
evolving  in  time,  it  would  convect  to  the  right,  with  the  vortex  lines  that  are  aligned  in  the 
axial  (z-direction)  direction  growing  longer.  There  would  also  be  self-induction  among 
all  the  lines.  For  now,  we  are  showing  the  basic  setup  and  ignoring  the  motion  of  the 
lines. 

These  vortex  lines  are  aligned  with  the  vortex  tube  in  the  axial  direction.  They 
emanate  from  the  vorticity  source.  However  at  the  end  of  the  vortex  tube,  they  branch  off 
radially  and  then  turn  back  upstream  where  they  end  at  the  z  =  0  boxmdary  that  defines 
the  location  of  the  vorticity  source.  It  can  be  helpful  to  think  of  these  radial  “spokes”  of 
vorticity  as  something  smiilar  to  the  starting  vortex  of  a  wing  s  bound  vortex  system. 

They  exist  to  insure  that  the  vortex  lines  do  not  end  in  the  fluid.  In  our  simulation,  we 
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will  allow  this  vortex  tube  end  and  these  radial  spokes  to  convect  far  enough 
downstream  that  they  will  have  a  negligible  effect  on  the  region  where  the  vortex 
breakdown  will  take  place.  Also,  we  will  assume  that  the  radial  spokes  extend  far 
enough  radially,  that  the  vortex  lines  that  have  turned  upstream  will  have  negligible  effect 
on  our  region  of  interest. 

However,  it  will  turn  out  that  a  gradient  in  vorticity,  brought  about  by  a  gradient 
in  the  z-direction  of  the  vortex  tube  circulation,  C,  will  be  important  to  our  study  of 
vortex  breakdown.  In  fact,  we  will  find  that  a  vorticity  gradient  will  initiate  the  vortex 
breakdown.  Thus,  we  must  create  this  circulation  gradient  without  violating  any  of  the 
Helmhotz  vortex  theorems. 

We  will  do  this  by  making  use  of  two  principles.  Figure  2-7  helps  to  describe 
these  ideas.  The  first  principle  is  the  superposition  of  vortex  lines.  We  can  linearly  add 
vorticity  fields  together,  just  as  we  add  induced  velocities  together  to  get  the  overall 
induced  velocity  on  each  vortex  filament  node.  Thus,  the  vorticity  source  will  feed  out 
additional  vortex  lines  at  z  =  0  at  the  times  we  specify  to  increase  the  circulation  in  the 
vortex  tube.  To  achieve  a  circulation  gradient  these  additional  vortex  lines  extend 
downstream  to  the  point  where  the  circulation  gradient  is  to  exist  and  then  branch  radially 
outward  in  radial  spokes  of  vorticity.  Just  as  we  used  the  radial  vorticity  spokes  to  bring 
the  vortex  lines  back  to  a  boimdary  where  they  met  the  end  of  the  vortex  tube,  we  will  do 
the  same  thing  to  insure  that  there  is  a  vorticity  gradient  in  the  z-direction. 
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Figxire  2-7  Changing  Circulation 
2.3. 1.1  Boundary  Conditions 

The  location  of  the  circulation  gradient  will  be  a  time-dependent  boundary 
condition,  which  we  will  explicitly  define  in  the  model.  Thus,  we  will  have  an  upstream 
boundary  condition  specified  at  the  vorticity  source.  There  will  be  no  boundary  condition 
downstream.  However,  we  will  allow  the  end  of  the  vortex  tube  to  convect  far  enough 
downstream  that  it  will  have  a  negligible  effect  on  the  vorticity  source  when  the 
circulation  gradient  is  introduced  into  the  flow.  Thus,  while  it  is  not  explicitly  specified, 
the  downstream  boundary  condition  is  a  vortex  tube  going  to  infinity.  This  is  inherently 
different  that  a  previous  study  of  vortex  breakdown  using  the  vortex  filament  method  by 
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Nakamura  et  al  (ref.  14),  where  they  specified  both  upstream  and  downstream 
boimdary  conditions  for  the  vortex  tube. 

Another  aspect  of  the  upstream  boundary  condition  is  that  we  will  only  introduce 
axial  vorticity.  However,  all  real  vortex  flows  typically  have  both  axial  and  azimuthal 
vorticity.  One  could  explicitly  specify  an  initial  azimuthal  vorticity,  as  do  Nakamura, 
Leonard,  and  Spalart  in  their  work  (ref  14).  In  doing  so  they  specify  an  initial  helical 
angle.  However,  we  wanted  to  allow  the  vortex  tube  to  assume  whatever  helical  angle  it 
naturally  would  (that  is,  whatever  helical  angle  the  induced  velocities  will  create).  Thus, 
we  specify  only  axial  vorticity,  and  allow  the  helical  angle  to  develop  on  its  own,  based 
on  the  azimuthal  velocity  induced  by  the  vortex  filaments  at  each  radial  location  in  the 
vortex  tube.  This  is  another  difference  from  the  Nakamura  et  al  study  where  they 
specified  the  helical  angle  of  the  vortex  tube.  Figure  2-8  illustrates  how  the  helical  angle 
can  develop  on  its  own.  Please  note  that  any  radial  vortex  lines  have  been  omitted  for 
clarity. 


Figure  2-8  Creating  Azimuthal  Vorticity 
The  filaments  are  in  a  fixed  azimuthal  and  radial  location  as  they  come  out  of  the 
vorticity  source.  The  specified  distribution  of  axial  vorticity,  ,  induces  an  azimuthal 

velocity.  Mg .  This  azimuthal  velocity  twists  the  vortex  lines,  creating  azimuthal  vorticity, 
cOg ,  by  means  of  vortex  lines  tilting. 

2.3.2  Numerical  case 

The  primary  numerical  simplification  is  to  not  include  vortex  filaments  that  are 
far  enough  away  that  their  induced  velocity  is  small  as  compared  with  the  induced 
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velocities  near  the  area  of  interest.  How  fer  is  far  enough?  There  are  two  ways  to 
determine  this.  The  first  way  is  to  look  at  two  cases.  One  case  will  have  the  far  away 
vortex  filaments  present  and  the  other  will  have  the  vortex  filaments  absent.  If  there  is 
little  or  no  change  between  the  two  cases,  we  can  be  sure  that  vortex  filaments  in 
question  can  be  safely  ignored. 

The  second  way  is  to  took  at  the  fiir  away  vortex  filaments  in  question  and  see  if 
they  are  affected  by  the  vortex  filaments  in  the  region  of  interest.  If  they  are  clearly  not, 
then  one  can  also  say  that  these  vortex  filaments  are  far  enough  away  to  safely  ignore 
them.  In  all  cases,  we  are  conservative  in  ignoring  vortex  filaments.  The  only  time  this 
comes  into  use  is  with  the  vortex  filaments  is  (a)  at  the  end  of  the  vortex  tube,  and  (b)  the 
vortex  filaments  in  the  radial  spokes  that  are  far  away  from  the  centerline  of  the  vortex 
tube.  For  (a),  in  our  numerical  model,  we  track  the  vortex  filaments  to  the  end  of  the 
vortex  tube,  which  we  have  allowed  to  travel  far  downstream.  However,  because  it  is  far 
downstream,  we  do  not  track  the  radial  spokes  and  upstream  traveling  vortex  lines  from 
the  end  of  the  vortex  tube.  For  (b),  we  track  the  radial  spokes  in  the  region  of  the 
vorticity  gradient  out  to  a  radius  of  approximately  ten  times  the  vortex  tube  radius.  At 
this  point,  we  track  them  no  further.  Nor  do  we  track  their  upstream  traveling 
components  of  vorticity. 

These  are  good  simplification  because  the  vortex  filaments  at  the  end  of  the 
vortex  tube  or  at  the  end  of  the  radial  spokes  from  the  circulation  increase  have  little 
influence  in  the  main  vortex  tube.  Hence,  we  can  infer  that  the  main  vortex  tube  is 
unaffected  by  them.  Since  it  would  make  no  difference  to  the  outcome  of  the  simulation 
whether  we  track  these  far  away  nodes  or  not,  for  numerical  efficiency,  we  choose  not  to 
track  them. 

Another  numerical  issue  is  that  of  how  we  model  the  vorticity  source  at  z  =  0.  It 
is  helpful  to  think  of  it  as  a  plane  where  vortex  filaments  with  a  specified  amount  of  axial 
vorticity  materialize.  Numerically,  we  treat  it  slightly  differently,  but  to  the  same  effect. 
In  the  numerical  model,  we  initially  prescribe  all  vortex  filaments  in  the  vortex  tube. 
However  we  locate  that  tube  entirely  upstream  of  the  z  =  0  plane.  The  z  =  0  plane 
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denotes  a  boundary  where,  to  the  left  (z  <  O),  which  is  upstream,  vortex  filaments  are 
only  moved  by  the  convection  of  the  fi-eestream  velocity,  .  Vortex  filaments  located 
at  z  <  0  are  not  allowed  to  induce  motion  to  each  other.  However,  as  soon  as  they  are 


convected  downstream  to  the  right  of  the  z  =  0  plane,  they  are  allowed  to  induce 
mutually,  in  addition  to  being  convected  downstream  by  the  free  stream  velocity.  This  is 
summarized  in  equations  (1-55)  and  (1-56)  and  illustrated  in  Figure  2-9. 


I 

z<0  z  =  0  z>0 

Figure  2-9  Numerical  construction  of  axial  vorticity  source 
A  further  numerical  simplification  we  make  is  to  exploit  symmetry  whenever  it  is 
present.  We  do  not  assume  the  flow  to  be  axisymmetric,  but  once  we  discover  that  it  is, 
we  can  employ  an  algorithm  that  allows  use  to  take  advantage  of  this  symmetry  to  reduce 
the  number  of  operations  and  decrease  model  run  time. 

2.4  Algorithm 

Now  let  us  look  at  the  general  algorithm  of  our  model. 

Figure  2-10  shows  an  outline  of  this.  We  will  look  at  each  part  specifically.  The 
complete  simulation  code,  written  in  C+h-  is  in  the  appendix. 
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Figure  2-10  Simulation  Algorithm 


2.4.1  Initialize  simulation 

The  first  step  is  to  initialize  the  simulatioa  We  specify  the  initial  location  of  all 
vortex  filament  segments.  We  assign  a  specified  amount  of  circulation  to  each  vortex 
filament.  Finally,  each  vortex  filament  segment  is  given  an  initial  radius  so  that  a  fixed 
vortex  filament  volume  can  be  calculated. 

2. 4. 1. 1  Specify  vortex  filament  location 

To  specify  the  location  of  the  vortex  filament  segments,  we  specify  the  location  of 
their  endpoint  nodes.  Each  node  is  the  endpoint  of  the  vortex  filament  segment,  so 
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specifying  the  node  location  defines  the  vortex  filament  segment  location.  In  the 
model  we  discretize  the  vortex  tube  into  nine  (9)  concentric  cylinders.  These  cylinders 
have  radii  of  0.42, 0.36, 0.3, 0.24, 0.18,  0.12, 0.08, 0.06,  and  0.03.  Each  cylinder  is 
further  discretized  into  sixteen  (16)  vortex  filaments,  which  are  evenly  distributed  around 
the  circumference  of  the  cylinder.  A  tracking  filament  is  also  specified  at  the  centerline, 
r  =  0.  These  centerline  vortex  filament  segments  contains  no  vorticity,  as  it  has  been 
distributed  to  the  surroimding  vortex  filaments.  Thus,  not  countmg  the  radial  spoke 
vortex  filaments,  there  are  9*16  +  1  =  1 45  vortex  filaments. 

Each  vortex  filament  is  then  discretized  axially  with  the  nodes  of  each  vortex 
filament  segment  initially  set  at  a  distance,  dz,  equal  to  0. 1 .  Specific  umts  are  not 
described,  but  they  are  units  of  length.  Enough  filaments  are  created  to  insure  that  there 
will  always  be  some  filaments  located  at  z  <  0  for  the  entire  length  of  the  simulation. 

Thus,  if  the  simulation  were  to  go  for  a  time  equal  to  1 .0  5  with  U„  =5  ,  the  vortex 

filaments  would  have  to  extend  back  to  at  least  z  =  -5.  This  would  require  at  least 

—  =  50  vortex  filament  segments  and  5 1  nodes  for  dz  =  0.1.  We  arrive  at  this  value 
dz 

through  an  optimization  process  that  will  be  discussed  in  a  later  chapter.  Thus,  each 
vortex  filament  node  has  its  own  location,  defined  by  f„  „  (x,  y,  z),  where  m  =  1  to  the 
maximum  number  of  axial  nodes  (iW)  and  n  =  1  to  the  maximum  number  of  vortex 
filaments  (N).  Figure  2-1 1  shows  how  r„„{x,y,z)  is  defined. 


Figure  2-11  Vortex  node  location 
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The  locations  are  initialized  according  to  equations  (2-3)  to  (2-6).  Though  we 
often  describe  things  in  cylindrical  coordinates,  Cartesian  coordinates  are  used 
throughout  the  actual  simulation. 

271 -(n-l)^ 


,„=r-cos|^- 
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(2-3) 
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(2-4) 


foT*  n  =j  to  N  =  16  •  y ,  for y  —  1  to  9,  for  m  —  I  to  M 
=z,init,  -dz-m 

j  ^2-5) 

for  n  =y  to  N  =  16  •  y ,  for y  =  1  to  9,  for  tw  =  1  to  M 

^m,n  ~  ’ y m,n (2“6) 


In  these  equations,  Tj  is  the  initial  radial  location  of  the  particular  concentric  of 

vortex  filaments,  and  z,initj  is  the  initial  location  of  the  end  of  the  vortex  tube.  When 

required,  radial  spoke  vortex  filament  segments  are  also  defined.  In  that  case,  the  nodes 
of  the  radial  spokes  are  given  the  z-location  where  the  circulation  increase  is  to  occur 
based  on  a  Galilean  transformation,  depending  on  what  time  the  circulation  increase  is  to 
come  out  of  the  vorticity  source.  Then,  they  are  assigned  radial  and  circumferential 
locations  to  coincide  with  the  existing  vortex  filaments  that  are  in  the  vortex  tube.  The 
radial  spoke  nodes  that  are  outside  of  the  vortex  tube  are  spaced  a  distance  dz  apart,  with 
whatever  number  of  vortex  filament  segments  required  to  place  the  ends  of  the  radial 
spoke  filaments  at  least  ten  times  the  radius  of  the  initial  vortex  tube.  In  this  case,  the 
vortex  tube  radius  is  0.42  and  the  radial  spoke  filaments  end  at  r  =  4.6.  The  following 
figures  (Figure  2-12  to  Figure  2-15)  show  location  of  all  the  vortex  filament  segments. 
For  clarity,  only  the  filaments  at  the  outermost  radius  and  the  radial  spoke  filaments  are 


shown. 
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Figure  2-12  Initial  location  of  vortex  lines,  perspective  view  one 
The  filaments  at  the  outermost  radius  of  vortex  filaments  in  the  vortex  tube  are 
shaded  to  look  like  a  solid  surface.  The  green  shaded  filaments  indicate  regions  of  the 
high,  final  value  of  circulation.  The  blue  regions  indicate  the  area  of  low  circulation. 
The  red  region  indicates  the  area  that  circulation  is  increasing.  Hence,  this  is  the  region 
where  the  radial  spokes  of  vorticity  exist  to  insure  solenoidality  of  the  vorticity  field. 
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Figure  2-13  Initial  location  of  vortex  lines,  perspective  view  two 
The  lines  that  run  along  the  length  of  the  vortex  tube  (in  the  axial  direction)  are 
the  actual  vortex  lines.  The  circumferential  lines  are  merely  tracking  lines  connecting 
nodes  of  vortex  segments  and  they  do  not  represent  vorticity  vectors.  In  this  initial 
configuration,  all  vortex  filaments  are  located  at  z  <  0 .  Thus,  they  are  not  allowed  to 
induce  each  other.  It  is  not  until  the  individual  filaments  are  convected  downstream  (+z) 
that  self-induction  will  be  allowed.  Figure  2-13  shows  a  closer  view  of  the  circulation 
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increase  region,  where  the  radial  spokes  extend  outward.  Once  the  radial  spokes  cross 
the  z  =  0  plane,  they  do  not  remain  purely  in  the  radial  direction  to  do  mutual  induction. 
This  region  is  sufficiently  upstream  of  the  end  of  the  vortex  tube  that  the  presence  of  that 
vortex  tube  end  has  no  bearing  on  the  induced  motion  in  the  circulation  increase  region. 


Streakiines,  Time  index:  0,  rinit:  0.42 


z,  axial-direction 

Figure  2-14  Initial  location  of  vortex  lines,  side  view 
Figure  2-14  shows  the  side  view  of  the  initial  vortex  filament  configuration.  In 
this  view  the  vortex  lines  run  from  left  to  right.  Interior  to  the  surface  visible  in  this 
figure  are  ninK  (9)  more  similar  cylinders,  each  composed  of  sixteen  (16)  axially  aligned 
vortex  filaments.  The  concentric  rings  are  clearly  visible  in  Figure  2-15,  which  shows  the 
initial  configuration  of  the  vortex  filaments  from  an  end  view.  Here  we  see  the  radial 
spokes  extending  outward.  We  also  see  nine  concentric  rings  of  vortex  filaments. 
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StreaWines,  Time  Index:  0,  rinrt:  0.03 


Figure  2-15  Initial  location  of  vortex  lines,  end  view 

2. 4. 1.2  Initialize  vortex  filament  core  radius 

Each  vortex  filament  has  the  cross  sectional  properties  of  a  real  vortex  tube.  Each 
vortex  filament  has  its  own  Gaussian  distribution  of  vorticity  and  its  own  core  radius. 
However,  because,  in  sum,  we  want  to  vortex  filaments  to  model  a  larger  vortex  tube,  we 
must  assign  each  vortex  filament  a  core  radius  that  creates  an  overlap  of  the  circulation  in 
each  vortex  filament.  This  overlap  will  give  the  entire  vortex  tube  a  smooth  vorticity 
distribution,  as  we  will  see  in  the  next  section.  Another  reason  we  assign  each  vortex 
tube  its  own  radius  is  that  it  allows  us  to  give  each  vortex  filament  segment  a  fixed 
volume.  By  continuity,  this  volume  will  remain  constant  throughout  the  simulation. 
However,  as  the  vortex  filament  segment  endpoints  are  stretched  apart  or  pushed  closer 
together  during  the  simulation,  the  core  radius  of  the  vortex  filament  segment  >vill  have  to 
change  to  maintain  the  constant  volume.  This  is  how  the  model  deals  with  the  stretching 
(or  contraction)  in  the  vortex  filaments.  If  the  nodes  are  stretched  apart,  the  core  radius 
of  the  vortex  filament  segment  will  have  to  decrease  to  maintain  the  constant  volume. 
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The  effect  in  the  simulation  is  that  the  velocity  induced  by  that  local  vortex  filament 
segment  increases  as  the  vorticity  increases.  The  circulation  remains  the  same,  but  it  is 
contained  in  a  smaller  cross  sectional  area.  Equation  (2-7)  shows  this  constraint. 
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m,n 


w+l,« 


-rr„,n\  =  V  =  const. 


(2-7) 


ntji 
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The  core  radius  of  the  vortex  filament  segment  on  the  ri^  vortex  filament  is 
a  .  The  vortex  filament  nodes  are  defined  by  f  j  and  „ .  The  value  of  is 

7n,n  *  ’ 

updated  at  each  time  step  according  to  Equation  (2-8).  In  our  model,  we  initially  set 
a  =  0.12  for  all  vortex  filaments.  This  radius  provides  sufficient  overlap  to  give  the 

fn,n 

overall  vortex  tube  a  very  smooth  distribution  of  vorticity.  Next,  we  will  assign 
circulation  to  each  of  these  vortex  filament  segments. 


2.4. 1.3  Assign  vorticity  to  each  filament 

After  the  location  of  the  vortex  filaments  has  been  specified  and  they  have  all 
been  given  a  core  radius,  vorticity  is  assigned  to  each  vortex  filament  segment.  The 
filament  segments  are  all  initially  aligned  in  the  axial  direction,  which  means  that  all  the 
vorticity  is  initially  aligned  in  the  axial  direction.  We  specify  vorticity  for  each  vortex 
filament  segment  by  giving  each  vortex  filament  segment  an  amount  of  circulation,  T,^^ . 

If  we  were  to  add  up  all  the  ^  at  a  particular  axial  location,  we  would  get  the  total 

circulation,  C,  at  that  location. 

Vorticity  will  remain  constant  along  a  vortex  line.  In  those  regions  where  the 
circulation  increases,  vortex  lines  of  constant  vorticity  are  superimposed  on  each  other, 
with  the  superimposed  vortex  lines  turning  radially  outward  at  the  appropriate  axial 
location  to  give  the  desired  total  circulation  profile  in  the  axial  direction.  However, 
vortex  filaments  at  different  radii  will  have  different  amounts  of  vorticity,  with  the 
vorticity  being  set  to  make  the  vorticity  distribution  interior  to  the  vortex  tube  match  that 
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of  an  ideal  vortex  tube.  Thus,  we  have  an  ideal  vorticity  distribution  that  we  want  to 
match.  In  this  case,  we  will  use  the  ideal  case  of  unbounded  swirling  flow.  For  this  case, 
the  circumferential  velocity  and  the  axial  vorticity  are  defined  by  Equations  (2-9)  and 
(2-10)  respectively. 


“e  = 


(2-9) 


\  d  (  \  C  b-e 

(d,=-—{ru,}  = - — 

r  dr  n 


-br^ 


(2-10) 


To  match  experimental  data,  we  choose  b  —  10.9  (ref.  14).  C  is  specified,  as  it  is 
the  total  circulation  at  a  given  axial  location.  We  cannot  exactly  match  a  continuous 
vorticity  distribution  in  a  discretized  vortex  tube,  but  we  can  do  it  in  a  way  that 


minimizes  the  error.  This  is  done  following  Nakamura,  Leonard,  and  Spalart  s  techmque 
using  the  Lagrange  multiplier  method  to  minimize  the  error  while  satisfying  the 
constraint  of  total  circulation  fixed  (ref  14).  Using  this  method,  we  achieve  the  initial 
axial  vorticity  profile  that  is  shown  in  Figure  2-16. 
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Figure  2-16  Ideal  vs.  actual  initial  axial  vorticity  profile 
As  was  stated  before,  this  total  circulation,  C,  at  a  given  axial  location  is  simply 
the  sum  of  the  circulation,  Fy  * . 

Figure  2-17  shows  the  total  circulation,  as  well  as  the  circulation  m  each 
concentric  ring  of  vortex  filaments.  To  say  it  another  way,  every  non-black  line 
represents  the  sum  of  the  circulation  from  the  sbcteen  vortex  filaments  in  each  of  the  nine 
rings  of  vortex  filaments.  The  sum  of  the  circulation  represented  by  these  non-black  lines 
in  the  total  circulation.  Thus,  we  see  peaks  of  circulation  near  the  radial  location  of  each 
ring  of  vortex  lines. 
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Figure  2-17  Circulation  sum 

The  axial  vorticity  in  each  concentric  ring  of  vortex  filaments  (this  is  the  vorticity  plotted 

in  each  non-black  line  in 


Figure  2-17)  is  given  by  equation  (2-11)  {ref.  14).  This  approximates  the 
distribution  of  axial  vorticity  in  each  vortex  filament  core  with  a  Gaussian  distribution. 


j  ~^'^j 


2ao]rXr] 


2  2  7 

:  "fr  H-aO) . 


for  y  =  1  to  9 


(2-11) 


In  equation  (2-1 1),  the  y-subscripted  values  refer  to  the  value  of  those  quantities  in 
an  entire  concentric  ring.  As  we  said  before,  the  vortex  tube  is  discretized  into  nine  (9) 
concentric  rings,  so  j  goes  fi’om  one  (1)  to  nine  (9).  C  is  the  circulation  of  the  entire 
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vortex  tube.  C  •  Fj  is  the  circulation  per  unit  circumference  of  the/*  ring  of  vortex 
filaments.  Thus,  C  has  units  of  circulation,  and  T,  has  units  of  [1/length].  We  choose  a 
=  0.413  to  achieve  a  Gaussian  distribution  of  vorticity.  The  circulation  in  an  individual 


271  •  r  •  C  •  r 

filament  in  that  ring  is - — - 


since  each  concentric  ring  is  made  up  of  sixteen 


(16)  individual  filaments.  Figure  2-18  shows  how  circulation  is  defined  for  a  single 
concentric  ring  of  vortex  filaments.  It  should  be  noted  that  the  filaments  overlap  both  in 
the  radial  and  in  circumferential  direction.  This  provides  the  continuous  vorticity 
distribution  that  is  seen  in  the  previous  figures. 


Figure  2-18  Circulation  in  a  ring  of  vortex  filaments 
Nine  (9)  rings  of  vortex  filaments  with  sbcteen  (16)  filaments  in  each  ring  are 
chosen  as  a  compromise  between  the  required  resolution  of  the  discretized  vortex  tube  to 
adequately  model  the  Gaussian  distribution  of  vorticity  in  the  vortex  core  and  numerical 
simulation  run  time.  The  following  table  shows  the  optimized  values  of  that 

271  rj  C  Tj  achieve  the  vorticity  distributions  in  the  preceding  figures. 
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In-r  -C- Tj 

Each  vortex  filament  at  radius  /*,  is  given  the  initial  circulation, - — - ,  as 

shown  in  this  table. 

Table  2-1  Vortex  filament  circulation  distribution _ 


r,  =  0.03 

r  ^  =  0.001640452531  •  C 

jMr,)  16 

rj  =  0.06 

r  ,,  ,=  ^  2  =0.002628594294  C 

jMh)  16 

rj  =  0.08 

r  ‘  ^  -  0.002429657471  •  C 

jMr,)  16 

r,  =  0.12 

r  =0.009090954102  C 

jMn) 

OO 

o 

11 

r  -0.01223904142- C 

r,  =  0.24 

r  ,,  ,  =  -0.01097097189-C 

jMr,)  ig 

=  0.30 

r  - 271  •  •  C  •  _  Q  01061059381  •  C 

Tg  =  0.36 

r  ^  ‘  -  0.004780605093  •  C 

rg  =  0.42 

2Tr-r  -C-T 

r  Z7t  L.  ^9  =0.008114549976-C 

7**^9  )  J  g 

2.4.2  Calculating  induced  velocities 

Now  that  the  vorticity  has  been  assigned  to  each  vortex  filament  segment,  we  can 
calculate  the  velocity  that  they  induce  on  each  other.  This  involves  a  decision,  though, 
because,  as  we  already  stated,  to  model  an  axial  vorticity  source,  we  will  only  allow 
vortex  filament  segment  endpoint  nodes  located  at  f„„{t)  =  to 

induce  velocities  or  to  receive  induced  velocities  if  they  are  located  at  z„  „  >  0 .  We  re¬ 
list  the  following  equations  to  show  the  distinction. 
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Z„,„<0  (2-12) 

^m,n  ^ induced, m,n (2-13) 


If  the  nodes  are  located  at  z  <  0 ,  then  they  only  receive  the  freestream  velocity, 
.  Likewise,  we  will  ignore  the  induced  velocities  of  any  vortex  filament  segment 


whose  endpoint  nodes  are  not  both  located  at  z  >  0 .  If  a  node  receives  only  ,  it  has 

all  the  velocity  that  it  will  get  during  the  time  step. 

Assuming  that  all  nodes  involved  in  the  following  equations  are  located  at  z  >  0 , 
we  will  proceed.  Recalling  the  discussion  of  the  Biot-Savart  law  in  the  previous  chapter, 
we  know  that  vorticity  at  point  P  will  induce  a  velocity  at  point  Q  according  to  equations 
(2-14)  and  (2-15).  Recall  that  the  s  function  allows  for  the  vortex  core  and  vortex 
stretching. 


(2-14) 


sw 


3 

2 


(2-15) 


Figure  2-19  Induced  velocity  from  P  on  Q 
In  our  model,  we  will  define  point  P  as  the  centroid  of  the  inducing  vortex 
filament  segment.  Point  Q  is  the  endpoint  node  of  a  vortex  filament  segment  that  is 
receiving  the  induced  velocity.  The  vortex  filament  core  radii,  c p  and  (Jq  ,  are  the 
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vortex  filament  radii  associated  with  the  inducing  vortex  filament  segment  at  point  P 
and  the  induced  vortex  filament  node  at  point  Q,  respectively.  Figure  2-20  shows  how 
we  apply  points  P  and  Q  to  vortex  filament  segments. 


Figure  2-20  Points  P  and  Q  on  vortex  filament  segments 
With  these  locations  defined,  we  can  calculate  the  induced  velocity  a  point  Q  due 
to  all  other  vortex  filaments.  We  do  this  according  equation  (2-16)  and  (2-17)  which  is 
simply  a  discretized  integration  of  the  Biot-Savart  law  fi’om  equation  (2-14). 


_  1  (r(0  X  A/(0)-  (0||,a^^„(0,a^.  t  (t)) 

~  .  2-i2^  M  ll:::,.l|3 


(2-16) 


for  w  =  1  to  J  for  n  =  1  to 


In  these  equations,  M„  „(Ois  the  velocity  induced  on  the  node  located  at 

vortex  filaments  at  time,  t.  The  limits  of 

summation  J  and  K  are  the  maximum  number  of  nodes  located  at  positions  z  >  0  and  the 
total  number  of  filaments,  respectively.  The  limit  is  J-1  because  there  are  J  nodes  and  J-\ 
vortex  filament  segments.  The  circulation,  T^^ ,  is  the  circulation  of  an  individual  vortex 

filament.  The  rest  of  the  terms  are  defined  in  equations  (2-18)  through  (2-27). 
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(2-18) 
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II 

(2-19) 

A/,(0=^(0;n.t-^(0y.u 

(2-20) 

^.(0  = 

^j+l,k  (0  ^  /+*•*  ~  ^m,n  (0 

(2-21) 

^,(0  =  1 

yj.ikiO +^{yj^u  (0  -  yj^i,k  (0)j  - 

(2-22) 

K(t)  = 

(0 + ^  (0  -  (0)j  -  (0 

(2-23) 

R(t)  =  [Rfit),R^{t),R,{t)) 

(2-24) 

(« (/)  X  a;(()),  =  A,  (0  •  A/,  -  (0  ■  A/, 

(2-25) 

(«(/)x4/(f)),  =«,(/)-Ai, 

(2-26) 

=«.(o-A/, 

(2-27) 

2.4.3  Update  vortex  filament  node  location 

With  the  induced  velocities  and  freestream  velocity  calculated  and  assigned  to 
each  node,  we  will  update  the  location  of  the  nodes  according  to  a  first  order  Euler 
integration,  as  in  equations  (2-28)  through  (2-30),  depending  on  the  location  of  .  In 

our  model,  we  choose  A/  =  0.01 .  We  arrive  at  this  value  through  an  optimization  process 
that  will  be  discussed  in  section  3.6 


for  z„,„<0; 

f„Jt  +  At)  =  r„^„it)+  j{UJ,)dt  ^„„it)  +  {Uj,yAt 

for  >0: 

+  ^0  =  r„Af)+  j{^ind.mAf)  +  UJ,)dt  ^„,„it)  +  (ui^,„,„(t)  +  UJ,yAt 

2.4.4  Update  Vortex  Filament  Segment  Core  Radius 

At  each  time  step,  the  core  radius  of  a  vortex  filament  segment  is  updated  to 
preserve  continuity  and  to  account  for  the  effects  of  stretching.  Recalling  that  we 
initialized  the  vortex  filament  segment  core  radius  according  to  the  following  equations, 
we  see  how  to  update  it. 


-  r„,M\\  =  V  =  const. 

(2-31) 

1  V 

(2-32) 
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(2-28) 

(2-29) 

(2-30) 


2.4.5  Repeat  main  program  loop 

With  the  new  node  position  and  the  updated  vortex  filament  segment  core  radii, 
we  can  update  the  time  step,  f  t  +  At ,  and  repeat  the  main  program  loop.  In  the  main 
program  loop,  we  output  the  current  vortex  node  location,  core  radii,  and  vorticity.  We 
then  calculate  the  velocity  on  each  vortex  node  from  every  vortex  filament  segment. 
Next,  we  update  the  position  and  core  radii  of  all  the  nodes  and  vortex  filament  segments. 
This  continues  until  t  =  t . 

2.5  Previous  vortex  filament  method  work 

Several  other  researchers  have  used  the  vortex  filament  method  to  study  vortex 
flows.  In  fact,  the  vortex  filament  method  is  really  just  an  extension  of  the  well  know 
vortex  lattice  method  for  calculating  aerodynamic  properties  over  three  dimensional 
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wings.  The  vortex  filament  method  takes  the  framework  of  the  vortex  lattice  method, 
adds  in  the  effects  of  a  viscous  vortex  core  and  vortex  stretching  and  tilting,  and  then 
removes  the  condition  that  the  vortex  must  remain  bound  to  the  wing.  In  the  vortex 
filament  method,  as  has  been  previously  described,  the  vortex  filaments  are  arranged  in  a 
vortex  tube,  allowed  to  induce  themselves,  and  then  advance  in  time  with  a  small  but 
finite  time  step. 

Two  previous  works  deserve  special  note  at  they  relate  to  this  research.  The  first  is 
Nakamura,  Leonard,  and  Spalart’s,  “Numerical  Simulation  of  Vortex  Breakdown  by  the 
Vortex-Filament  Method”  (ref.  14).  The  computational  method  for  our  present  study  is 
largely  based  on  their  work.  Specifically,  we  use  their  approximation  for  the  viscous  core 
and  their  method  of  taking  into  the  account  of  stretching.  Furthermore  we  base  the 
method  for  calculating  and  applying  the  initial  distribution  of  vorticity  to  match  the  cross- 
sectional  Gaussian  vorticity  distribution. 

Our  work  differs  in  model  resolution  and  the  boundary  conditions.  As  a  result,  the 
outcome  of  our  simulation  is  significantly  different.  In  terms  of  model  distribution,  due 
to  increases  in  computer  capabilities,  our  vortex  tube  is  discretized  into  many  more 
vortex  filament  segments.  In  terms  of  the  boundary  conditions,  Nakamura,  Leonard,  and 
Spalart  fix  both  the  upstream  and  downstream  ends  of  the  vortex  tube.  Additionally,  they 
specify  an  initial  helical  angle.  Furthermore,  they  give  the  vortex  tube  a  vorticity 
distribution  that  is  constant  along  the  length  of  the  vortex  tube.  Our  research  models  an 
axial  vorticity  source  where  the  vortex  tube  is  fed  downstream  by  a  freestream  velocity. 
The  amount  of  vorticity  created  at  the  vorticity  source  can  change  in  time  using  radial 
spokes  of  vorticity.  Figure  2-21  shows  the  results  of  the  Nakamura,  Leonard,  and  Spalart 
simulation  as  compared  to  the  present  results. 


Streaklines,  Time  index:  180,  Final  C:1S 


Figure  2-21  Nakamura,  Leonard,  and  Spalart  results  (top),  new  results  (bottom) 

The  second  work  we  should  note  is  that  of  Saghbini  and  Ghoniem  in  their  paper 
“Numerical  Simulation  of  the  Dynamics  and  Mixing  in  a  Swirling  Flow”  {ref.  15).  As 
we  were  formulating  our  ideas  we  came  across  their  work.  They  use  a  sunilar  vortex 
filament  method  to  study  vortex  breakdown.  Much  like  our  model,  they  simulate  a 
vorticity  source  that  convects  downstream,  without  an  explicit  downstream  boundary 
condition.  However,  several  things  are  not  clear.  The  first  is  how  they  can  account  for 
the  effects  of  diffusion  while  still  using  a  vortex  line  method  without  violating  the 
Helmholtz  vortex  theorem  that  vorticity  follows  the  fluid  element.  Others  have  said  that 
this  limitation  is  “one  of  the  biggest  obstacles  toward  extending  it  [the  vortex  filament 
method]  to  viscous  flows”  {ref.  16).  Furthermore,  they  apparently  use  the  finite  length  of 
the  vortex  tube  as  the  source  of  a  vorticity  gradient  that  leads  to  vortex  breakdown.  They 
state  that  “[a]t  time  zero,  a  negative  axial  vorticity  gradient  is  established  due  to  the  finite 
length  of  the  vortex  lines”  {ref.  15).  While  this  vorticity  gradient  does  seems  to  be  a 
necessary  condition  for  vortex  breakdown,  it  seems  improper  to  model  it  in  a  way  that 
violates  another  one  of  Helmhotz’s  vortex  theorems,  namely  that  a  vortex  line  cannot  end 
in  a  fluid.  Numerically,  one  can  say  that  the  end  of  a  finite  length  vortex  line  is  far 
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enough  away  from  the  region  of  interest  that  one  can  say  that,  for  all  practical 
purposes,  the  vortex  line  may  as  well  extend  downstream  to  infinity.  It  is  not  clear  that 
this  is  the  case  here.  Furthermore,  several  asymmetries  appear  to  develop  in  their 
solution  that  are  not  specifically  accounted  for.  In  a  vortex  filament  model,  if  everything 
is  initially  symmetric  in  the  flow,  and  the  algorithm  is  applied  properly,  it  should  always 
remain  symmetric.  Whether  things  were  initially  symmetric  or  not  in  their  model  is  not 
clear. 

Still,  Saghbini  and  Ghoniem  bring  forward  the  idea  about  the  importance  of  the 
vorticity  gradient.  They  point  to  the  vorticity  gradient,  coupled  with  radial  expansion  as 
the  initiating  cause  of  vortex  breakdown.  Furthermore,  they  identify  the  azimuthal 
vorticity  sign  switch,  first  supposed  by  Brown  and  Lopez  {ref.  17),  as  being  part  of  a 
feedback  mechanism  that  leads  to  a  steady  state  vortex  breakdown.  We  largely  agree, 
and  we  reach  similar  conclusions  with  our  model,  thus  adding  to  the  evidence  that 
supports  the  validity  of  these  claims.  We  will  attempt  to  be  more  rigorous  with  the 
Helmholtz  vortex  theorems. 
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3.  Vortex  Filament  Method  Results 

3.1  Model  Output 

Before  we  present  out  test  cases,  we  need  to  explain  exactly  what  the  oxitput  of 
those  cases  shows.  First,  we  need  to  discuss  the  differences  and  similarities  between 
streaklines,  pathlines  as  they  relate  to  materials  lines,  and  vortex  lines  in  unsteady  flows. 
Then  we  will  present  an  example  of  the  results  and  show  what  that  output  describes. 

3.1.1  Streaklines  and  pathlines 

In  flow  visualization  experiments,  the  researcher  will  often  uses  tracers,  be  they 
smoke  or  dye,  to  view  the  path  that  the  particles  take  in  whatever  flow  is  being  studied. 
The  lines  that  these  tracers  form  are  called  streaklines.  We  call  them  streaklines  because 
the  tracers  often  do  not  appear  as  distinct  particles,  but  in  a  continuous  streaking  line.  If 
we  were  to  mark  the  vortex  over  a  delta  wing  in  a  water  tunnel  with  dye,  we  would  see 
the  picture  that  is  in  Figure  3-1 


Figure  3-1  Streaklines  in  a  vortex  flow 
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Though  we  are  marking  particles,  we  cannot  distinguish  the  actual  movement 
of  a  particular  particle.  In  steady  flow,  this  presents  no  difficulty  because  the  pathline 
and  the  streakline  are  coincident.  However,  in  periodic  regions  of  steady  flow,  or  in  all 
regions  of  unsteady  flow,  pathlines  and  streaklines  are  not  coincident  (ref.  18).  This  can 
be  easily  illustrated  in  the  following  example,  illustrated  in  Figure  3-2. 


Figure  3-2.  Pathlines  vs.  streaklines 


Imagine  we  have  a  chimney  that  sends  out  smoke.  There  is  breeze  that  blows  the 
smoke  to  the  right.  At  time  Xi  the  wind  changes  direction.  Here  we  can  clearly  see  the 
difference  between  the  actual  pathlines  of  the  particles  and  the  streaklines.  Streaklines 
show  all  particles,  released  successively  from  a  fixed  point,  at  one  instant  in  time:  Po,  Pi, 
P2,  P3,  at  time  t4.  Pathlines  shows  how  one  particle  travels  over  time:  Po  at  to,  ti,  t2,  ts,  and 

t4- 

In  vortex  flows  the  vortex  breakdown  is  inherently  periodic  or  unsteady,  and  the 
fact  that  the  streaklines  we  see  in  flow  visualization  experiments  do  not  reveal  the 
pathlines  has  caused  much  difficulty.  One  thing  that  can  happen  is  that  a  region  of  flow 
may  appear  to  rotate  in  the  opposite  direction  of  the  rest  of  the  flow,  when  in  reality  it  is 
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rotating  in  the  same  direction,  only  more  slowly.  Particles  that  are  actually  moving  on 
a  very  benign  path  may  appear,  by  looking  only  at  their  streaklines,  to  move  violently  and 
in  odd  directions.  This  illusion  is  especially  strong  in  the  region  downstream  of  vortex 
breakdown.  Therefore,  one  must  be  very  careful  when  looking  at  streaklines  not  to  draw 
too  many  inferences  about  the  actual  fluid  motion.  To  do  this  we  need  to  look  at  velocity 
and  vorticity  fields  or  we  need  to  track  individual  particles.  With  the  vortex  filament 
method,  velocity  and  vorticity  fields  are  associated  with  the  location  of  all  the  vortex 
filaments  at  a  given  instant  in  time.  These  vortex  filaments  can  be  regarded  as  streaklines 
in  that  they  show  the  position  of  many  vortex  filaments  segments  or  fluid  particles, 
‘released’  at  a  fix  position,  z  =  0,  at  a  single  instant  in  time.  Thus,  we  can  use  their 
position  to  discuss  induced  velocity  and  associated  vorticity,  but  we  can  not  use  their 
relative  positions  at  any  one  instant  in  time  to  say  much  of  anything  about  the  fluid 
motion  over  time. 

5. 1.2  Vortex  filament  method fiow  visualization 

At  this  point,  we  will  not  discuss  what  any  output  tells  us  about  the  nature  of 
vortex  breakdown.  We  will  only  discuss  how  we  will  present  that  output.  Remember 
that,  in  this  model,  we  discretized  the  vortex  tube  into  nine  (9)  concentric  rings  of 
vorticity.  We  fiirther  discretized  each  ring  into  sixteen  (16)  vortex  filaments  plus  an 
additional  center  tracking  line.  As  we  previously  discussed,  the  vortex  lines  are 
equivalent  to  streaklines.  Thus,  what  we  see  in  the  model  output  corresponds  to  dye 
streaks  in  experimental  flow  visualization,  as  in  Figure  3-1.  However,  if  we  show  all  the 
vortex  filaments  segments,  we  will  see  a  very  tangled  picture  of  vortex  filaments.  Figure 
3-3  and  Figure  3-4  show  a  perspective  view  of  the  output  with  all  the  vortex  filaments  at 
T  =  1.8.  Again,  we  will  discuss  the  conditions  the  generated  this  case  later.  For  now  we 
are  only  concerned  with  the  presentation  of  the  output.  The  fi'eestream  velocity  always 
acts  in  the  positive  axial  direction.  We  will  indicate  this  on  the  next  two  figures  with 
arrows,  and  assume  it  fi'om  then  on  with  the  arrows  omitted. 
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Streaklines,  Time  index:  180,  Final  C:15 


Figure  3-3  Example  of  model  output,  all  vortex  filaments,  side  view 


Streaklines,  Time  index:  180,  Final  C;15 


Figure  3-4  Example  of  model  output,  all  vortex  filaments,  perspective  view 
As  you  can  see,  it  is  difficult  to  see  the  structure  of  the  vortex  tube  when  it  is 
described  by  a  mass  of  lines.  We  have  even  left  out  the  radial  spokes  of  vorticity  (they 
are  present  in  the  simulation  in  all  cases,  but  we  may  leave  them  out  of  the  output  for 
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clarity).  Figure  3-5  shows  the  same  perspective  view  with  the  radial  vortex  spokes 
included.  As  you  can  see,  the  results  get  even  more  muddled. 


Figure  3-5  Example  of  model  output,  with  radial  spokes,  perspective  view 
To  make  the  output  clearer,  we  will  often  only  show  one  or  two  of  the  discretized 
vortex  tube  radii.  We  will  often  include  the  center  tracking  line  because  it  gives  clear 
picture  of  the  relative  flow  at  the  very  center  of  the  vortex  tube.  Thus,  in  the  following 
figures,  we  include  the  r  =  0.42  and  the  r  =  0.12  concentric  sleeves  of  vortex  filaments. 
Figure  3-6  shows  the  side  view  and  Figure  3-7  shows  a  perspective  view. 
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Streaklines,  Time  index:  180,  Final  C:15 


Figure  3-6  Example  of  model  output,  r  =  0.42,  r  =  0.12,  side  view 
Streaklines,  Time  index:  180,  Final  C:15 


9  -1 


Figure  3-7  Example  of  model  output,  r  =  0.42,  r  =  0. 12,  perspective  view 
By  just  removing  all  but  two  radii,  the  picture  becomes  much  clearer.  However, 
even  with  two  radii  it  can  be  difficult  to  see  specific  features  of  the  output.  So  we  may 
want  to  look  at  just  one  radius,  as  in  Figure  3-8  and  Figure  3-9. 
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Streaklines,  Time  index:  180,  Final  C:15 


Figure  3-8  Example  of  model  output,  r  =  0.12,  side  view 
Streaklines,  Time  Index:  180,  Final  C;15 
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Figure  3-9  Exanq)le  of  model  output,  r  =  0.12,  perspective  view 


This  is  the  output  that  most  closely  matches  what  one  would  see  in  a  dye  streak 
flow  visualization  experiment.  There  dye  gets  entrained  into  a  single  interior  radial  layer 
of  the  vortex  tube  and  remains  with  that  layer,  marking  its  vortex  lines  to  form 
streaklines.  Experimentally,  you  can  mark  different  radial  layers  of  the  vortex  tube,  or 
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you  can  just  mark  one.  The  previous  figures,  where  we  show  only  one  ring  of 
concentric  vortex  filaments,  are  equivalent  to  that.  All  the  other  vortex  filaments  are 
present  in  the  simulation,  just  as  they  are  present  in  a  real  vortex  tube.  They  affect  the 
movement  of  the  marked  radius,  but  we  only  see  that  one  marked  radius,  just  as  in  the 
one  radius  of  vortex  filaments  that  we  are  showing. 

Even  this  view,  while  better  than  the  others,  has  its  difficulties.  For  instance,  in 
Figure  3-8,  the  side  view,  it  is  not  clear  which  way  the  vortex  filaments  are  tilting 
because  we  seen  vortex  lines  on  both  sides  of  the  vortex  tube.  To  overcome  this,  we  will 
render  the  vortex  lines  as  a  surface,  with  the  axial  running  lines  being  the  vortex 
filaments  and  the  circumferential  lines  being  simply  lines  that  mark  the  axial  length  of  a 
vortex  filament  segment.  Figure  3-10  and  Figure  3-11  give  a  side  and  perspective  view 
of  this  rendered  case. 


Streaklines,  Time  index:  180,  rinit:  0.12 
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Figure  3-10  Example  of  model  output,  r  =  0.12,  rendered  side  view 
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Streaklines,  Time  index:  180,  rinit;  0.12 


Figure  3-1 1  Example  of  model  output,  r  =  0.12,  rendered  perspective  view 
We  can  use  a  closer  view  to  look  at  the  specific  features  of  this  rendered  view. 
Figure  3-12  shows  one  such  image.  In  this  view,  the  yellow  line  shows  one  particular 
vortex  filament.  The  purple  line  shows  the  circumferential  spacing  line  that  indicates  the 
length  of  the  vortex  filament  segments.  Each  intersection  of  a  vortex  filament  segment 
and  a  circumferential  tracking  line  indicates  a  vortex  filament  segment  node.  The  green 
shaded  regions  show  areas  of  the  vortex  tube  with  a  large  circulation.  The  blue  shaded 
regions  show  areas  of  the  vortex  tube  with  low  circulation,  and  the  red  shaded  regions 
show  areas  where  the  circulation  of  the  vortex  tube  is  changing  due  to  the  superposition 
of  vortex  lines.  We  will  describe  why  we  change  the  circulation  along  the  vortex  tube  in 
the  next  section. 
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Streaklines,  Time  index:  180,  rinit:  0.12 


4.5  -0.5 


Figure  3-12  Closer  view,  r  =  0.12,  rendered  perspective  view 
An  intrinsic  benefit  to  the  vortex  filament  method  is  to  be  able  to  infer  induced 
velocities  by  looking  at  the  orientation  of  the  vortex  filaments.  We  know  that  the 
vorticity  vector  is  always  aligned  with  the  vortex  filaments,  so  simply  by  using  the  right 
hand  rule,  we  can  quickly  see  what  the  induced  velocity  should  be  from  any  given  vortex 
filament.  This  provides  insight  into  how  vortex  breakdown  occurs  that  cannot  be  found 
by  looking  more  complicated  vortex  simulations  techniques.  These  techniques  do  not 
intrinsically  preserve  the  vortex  line,  so  often  plots  of  the  vortex  filaments  are  not 
present.  To  see  the  loss  in  information,  we  can  compare  the  vortex  filament  simulation 
results  with  and  without  the  vortex  lines  present.  Figure  3-13  shows  the  rendered  results 
at  T  =  1 .8.  While  we  can  see  the  vortex  breakdown  in  the  top  image,  we  can  see  the 
vortex  filaments  and  determine  the  resulting  induced  velocities  that  lead  to  this 
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breakdown  in  the  bottom  image.  The  ability  to  visualizing  these  vortex  filaments  is 
one  of  the  principal  benefit  of  the  vortex  filament  method. 

Vortex  Filaments,  Time:  1 .8 


4.5 


Figure  3-13  Simulation  results,  with  (bottom)  and  without  (top)  vortex  filaments,  T  -  1.8 
3. 1.3  Vorticity  and  Velocity  Data 

Besides  presenting  the  vortex  filament  locations  at  a  given  instant  in  time,  we 
have  also  calculated  the  velocity  and  vorticity  associated  with  the  vortex  tube  on  the 
meridional  plane  passing  through  the  vortex  tube.  This  was  done  by  simply  taking  the 
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vortex  filaments  at  a  given  point  in  time,  with  their  circulation  and  vortex  core  radii 
known,  and  calculating  their  induced  velocity  on  fixed  nodes  evenly  placed  on  the 
meridional  plane  of  the  vortex  tube.  Figure  3-14  and  Figure  3-15  show  illustrations  of 
this.  The  first  shows  a  perspective  view  with  the  induced  azimuthal  velocity  presented  on 
the  meridional  plane.  The  second  shows  a  side  view  of  the  vortex  tube  overlaid  with  the 
axial  velocity  (including  the  fi-eestream  velocity)  on  the  meridional  plane. 

Circumferentiai  Velocity,  meridtona)  plane,  Time  index:  1 80 


9  -15 

Figure  3-14  Example  of  azimuthal  velocity,  meridional  plane,  perspective  view 
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Figure  3-15  Example  of  axial  velocity,  meridional  plane 


The  velocity  is  in  units  of  [lengfh\l[tme]  and  the  vorticity  is  in  units  of  \l[time\ 
Specific  units  are  not  used  because  we  are  concerned  with  the  differences  in  values  and 
not  their  absolute  value.  From  this  velocity  field,  we  can  calculate  the  vorticity  field. 
Calculating  the  axial  and  radial  vorticity  requires  that  we  calculate  induced  velocities 
slightly  away  fi'om  each  side  of  the  meridional  plane  in  addition  to  the  meridional  plane 
velocities.  We  calculate  the  vorticity  according  to  the  second  order  scheme  in  equations 
(3-l)and  (3-2)  (r^  19).  This  equation  shows  the  calculation  of  co^.  The  calculations  for 


coe  and  cOr  are  equivalent.  In  the  meridional  plane,  eoe  =  and  =  (Ox- 


■  '.y 


4Ax- Ay 


(3-1) 


^(^x)i,j+t  (^x)i-IJ+l) 

The  equations  are  applied  to  evenly  distributed  nodes  on  (or  slightly  ofi^  the 
meridional  plane  where  we  have  calculated  the  induced  velocity. 


(3-2) 
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3.2  Application  to  vortex  breakdown 

Much  of  the  current  theory  on  vortex  breakdown  deals  with  the  breakdown  in  a 
steady  state  condition.  This  theory  is  very  useful.  However,  it  must  typically  assume  an 
initiating  condition  in  the  vortex  tube,  such  as  a  radially  expanded  vortex  tube  {ref.  17). 
Often,  a  pressure  gradient  is  assumed  to  initiate  the  breakdown  {ref.  20).  We  want  to  take 
a  step  back  fi-om  that  and  look  at  the  unsteady  formation  of  vortex  breakdown,  starting 
with  a  constant  radius  vortex  tube  fed  fi'om  a  vorticity  source  and  seeing  what  input 
conditions  are  required  to  initiate  breakdown.  Vortex  breakdown  may  not  require  a 
pressme  gradient,  but  it  may  be  brought  about  through  a  self-induction  caused  by  a 
vorticity  gradient  in  the  vortex  core.  Initial  research  into  a  possible  self-induction  method 
that  initiates  vortex  breakdown  {ref.  21  and  22),  supports  this  conclusion  for  the 
importance  of  a  vorticity  gradient.  We  conclude  that  a  vorticity  gradient,  with  axial  and 
azimuthal  vorticity  in  the  vortex  tube  increasing  in  the  upstream  direction,  is  an  important 
initiator  of  the  vortex  breakdown  mechanism.  Others  have  concluded  this  as  well  {ref. 
15).  Therefore,  we  will  use  the  previously  described  vortex  filament  method  to  test  this 
hypothesis.  If  a  vorticity  gradient  leads  to  a  configuration  of  the  vortex  tube  that  has  the 
same  conditions  as  a  vortex  breakdown  in  a  real  flow,  we  can  conclude  that  a  vorticity 
gradient  can  initiate  vortex  breakdown.  Furthermore,  with  the  vortex  filament  method, 
we  can  observe  how  the  vorticity  gradient  changes  the  vortex  tube  so  that  we  might 
understand  the  specific  imsteady  mechanism  that  leads  to  a  steady  vortex  breakdown. 
With  this  insight,  one  can  devise  specific  ways  to  prevent,  weaken,  or  delay  vortex 
breakdown  in  real  flows. 

In  our  model,  we  will  generate  a  vorticity  gradient  by  causing  a  change  in 
circulation  in  the  vortex  tube.  This  is  not  the  only  way  a  vorticity  gradient  can  be  created 
in  a  vortex  tube.  It  is  possible  that  a  pressure  gradient  causes  a  vorticity  gradient  in  the 
vortex  tube,  which  would  then  initiate  vortex  breakdown.  This  may  be  a  fine  point,  but  it 
shows  that  the  physical  mechanism  for  vortex  breakdown  is  present  in  the  vortex  tube 
itself.  After  we  present  the  results,  we  will  discuss  other  possible  ways  to  achieve  this 
vorticity  gradient,  especially  in  vortex  tube  over  real  wings. 
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3.3  Vorticity  gradient  in  the  discretized  vortex  tube 

In  this  case,  we  will  cause  a  vorticity  gradient  in  the  vortex  tube  by  specifying  a 
change  in  circulation  in  that  vortex  tube.  As  was  previously  described,  we  change  the 
circulation  in  the  vortex  tube  by  the  superposition  of  vortex  filaments  that  turn  outward  in 
radial  spokes  of  vorticity  at  the  required  axial  location  to  achieve  the  prescribed 
circulation  profile  while  maintaining  solenoidality  within  the  vortex  tube.  For  a  negative 
vorticity  gradient  case  (vorticity  in  the  vortex  tube  increasing  upstream),  the  circulation 
in  the  vortex  tube  coming  out  of  the  vorticity  source  as  a  fimction  of  time  is  given  the 
profile  in  Figure  3-16.  The  circulation  has  the  unspecified  units  oi[lengthfl[time\  This 
is  the  circulation  of  the  vortex  tube  as  in  leaves  the  vorticity  source.  Those  filaments 
leaving  at  that  time  keep  their  given  circulation  for  the  entire  simulation.  Therefore, 
downstream  of  the  vorticity  source,  there  may  be  vortex  filaments  with  a  lower 
circulation  than  those  leaving  the  vorticity  source  at  the  current  instant  in  time. 


0  0.2  0.4  0.6  0.8  1  1.2  1.4  1.6  1.8 

T,  time 

Figure  3-16  Circulation  of  vortex  tube  at  vorticity  source 
For  the  circulation  profile  in  Figure  3-16,  the  initial  circulation  of  the  vortex  tube 
is  very  low  (C  =  0.1  [lengthf  l[time^.  It  then  starts  to  increase  (as  it  is  fed  out  of  the 
vorticity  source)  at  T  =  0.6  [time\  The  vorticity  source  continues  increasing  its 
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circulation  output  until  T  =  0.8  [time]  when  the  circulation  reaches  its  final  value  (C  = 

15  [lengthfl[tme]).  In  this  case,  U^=5  [length]/[time].  The  end  of  the  vortex  tube  is 

far  enough  downstream  not  to  affect  the  region  of  interest. 

One  should  remember  that,  while  we  talk  about  the  vorticity  source  changing  its 
circulation  output  in  time,  the  vortex  tube  has,  in  the  actual  model  algorithm,  already 
been  assigned  the  circulation  increase  in  the  region  of  z  <  0 .  However,  it  is  only  able  to 
affect  the  vortex  tube  when  the  vortex  filaments  with  the  higher  circulation  convect 
downstream  to  a  location  z  >  0 .  Therefore,  we  can  speak  either  in  terms  of  the  vorticity 
source  changing  its  circulation  output  in  time,  or  in  terms  of  the  vortex  filaments  pre¬ 
assigned  in  z  <  0  convecting  into  the  z  >  0  region.  It  is  more  physical  to  talk  about  the 
former  case,  so  that  is  how  we  will  describe  it. 

3.4  Primary  Case 

We  will  present  the  case  of  the  simulation  responding  to  the  circulation  change 
described  in  Figure  3-16.  We  will  present  the  vortex  filament  plots,  as  previously 
described,  in  a  variety  of  forms  to  convey  the  results.  We  will  present  velocity  and 
vorticity  color  contour  plots  along  the  meridional  plane  for  the  entire  vortex  tube.  The 
following  chart  summarizes  the  case  parameters.  We  will  present  the  results  at  intervals 
ofT  =  0.3. 


Table  3-1  Primary  case  parameters 


Az,  initial  vortex  filament 

segment  length 

0.1  [lengfh] 

At,  simulation  time  step 

0.01  [time] 

c 

^  initial 

0.1  [lengtW^I[time] 

^ final 

15  [lengthfl[time] 

Tc 

'•'initial 

0.6  [time] 

Tc 

'•'final 

0.8  [time] 
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3. 4. 1  Vortex  filament  plots 

We  present  four  different  plots  of  the  simulation  results.  The  first  shows  a  side 
view  of  the  vortex  filaments  that  start  at  r  =  0.12.  The  second  set  of  plots  show  the  same 
view  with  both  r  =  0.12  and  r  =  0.42  vortex  filaments.  The  third  set  includes  both  r  = 
0.12  and  r  =  0.42  and  adds  the  radial  spokes.  The  last  set  shows  a  side  rendered  view  of 
the  r  =  0.12  vortex  filaments. 


Filaments,  Time:  2.4 


Vortex  filaments,  r  =  0.12,  r  =  0.42 

Vortex  Filaments,  Time:  0 
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Vortex  Filaments,  Time:  0.6 
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Vortex  RIaments,  Time:  1 .8 


5  6 

Filaments,  Time:  0.6 


5  6 

Filaments,  Time:  1.2 
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5. 4.2  Velocity  and  Vorticity  plots  on  meridional  plane 

In  the  following  plots,  we  present  the  axial  and  azimuthal  velocity  and  vorticity, 
as  calculated  along  the  meridional  plane  of  the  vortex  tube.  The  vortex  filaments  at  r  = 
0.12  and  r  =  0.42  are  overlaid  to  show  how  their  instantaneous  location  correlates  with 
the  vorticity  and  velocity  data.  A  color  bar  scale  is  located  at  the  end  of  each  section. 

3. 4. 2. 1  Color  contour  axial  velocity  plots 

Axial  Velocity,  meridonal  plane,  Time:  0 
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Axial  Velocity,  meridonal  plane,  Time:  0.6 


3. 4. 2. 2  Color  contour  azimuthal  velocity  plots 


Azimuthal  Velocity,  meridonal  plane,  Time:  0 


Azimuthal  Velocity,  meridonal  plane,  Time:  0.6 


Azimuthal  Velocity,  meridonal  plane,  Time:  1.2 
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3. 4.2.4  Color  contour  axial  vorticity plots 


Axial  Vorticity,  meridional  plane,  Time:  0 
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Axial  Vorticity,  meridional  plane,  Time:  1 .8 
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3J.3  Streamline  plot 

Based  on  the  meridional  velocity  data,  we  can  construct  a  plot  of  the  streamline  at 
the  meridional  plane.  The  plot  below  shows  the  streamlines  at  T  =  2.4.  The  main 
recirculation  region  is  clear,  but  another  area  of  radial  expansion  is  beginning  to  form 
downstream  of  the  main  breakdown.  While  not  as  clear  in  the  vortex  filament  plots,  this 
secondary  breakdown  region  is  very  clear  in  this  streamline  plot. 

2 
1.5 
1 

0.5 
0 

-0.5 
-1 
-1.5 
-2 

Figure  3-17  Streamlines,  T  =  2.4,  meridional  plane 
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3.4.4  Case  1  pathlines,  rendered 

As  was  previously  stated,  with  streaklines,  very  odd  patterns  can  occur  in  the  flow 
visualization,  while  the  actual  particle  pathlines  may  appear  very  different.  The 
following  plots  show  the  pathlines  for  Case  1,  for  particles  released  at  the  below  stated 
radii  and  the  stated  time.  In  all  cases,  the  final  time  is  T  =  2.7. 

3.4.4. 1  Pathlines,  ri„it  =  0.12 


Pathlines,  Initial  Radius:  0.12,  Particles  Released  at  Time;  0 


Pathlines,  Initial  Radius:  0.12,  Particles  Released  at  Time;  0.6 
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Pathlines,  Initial  Radius;  0.12,  Particles  Released  at  Time:  1.5 
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3. 4. 4.2  Pathlines,  ri„u  =  0.42 

Pathlines,  Initial  Radius;  0.42,  Particles  Released  at  Time:  0 


1.5f 
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Pathlines,  Initial  Radius;  0.42,  Particles  Released  at  Time:  1 .5 


101 


3.5  Other  Cases,  Vary  Cfmai 

The  following  results  show  how  the  output  changes  when  the  final  circulation  (see 
Figure  3-16),  and  hence  the  strength  of  the  vorticity  gradient,  varies  from  C  =  1  to  C  = 

25.  The  results  are  presented  at  T  =  1.8.  The  initial  circulation  remain  C  =  0.1.  It  is 
clear  from  these  results,  that  the  higher  the  final  circulation,  the  sooner  the  characteristics 
of  vortex  breakdown  develop.  This  points  to  some  type  of  criteria  involving  the  strength 
of  the  vorticity  gradient  whether  vortex  breakdown  will  develop  over  a  wing  or  not.  If 
the  critical  factor  in  vortex  breakdown  is  the  formation  of  the  turning  point  (the  turning 
point  is  where  the  vortex  filaments  turn  inward  on  themselves,  see  section  4.4.4  for 
details),  then  we  can  say  the  following.  If,  for  a  given  negative  vorticity  gradient  in  a 
vortex  tube,  the  turning  point  forms  at  T  =  Ttp,  then  if  Ttp  is  greater  than  the  time  for  the 
region  of  that  vorticity  gradient  to  convect  downstream  over  the  surface  of  the  wing,  then 
a  steady,  sustained  vortex  tube  will  not  develop.  Conversely,  if  Ttp  is  less  than  the  time  it 
takes  for  the  region  of  the  vorticity  gradient  to  convect  downstream  of  the  wing  surface, 
there  will  be  a  steady,  sustained  vortex  breakdown  over  that  vwng.  Other  factors  such  as 
a  velocity  gradient  and  the  vortex  sheet  that  constrains  the  vortex  core  would  also  need  to 

be  taken  into  account  to  eloborate  fuller  on  this  idea. 

Vortex  Filaments,  C, final:  1 ,  Time:  1 .8 


-0.5 


3  4  5 

Vortex  Filaments,  C, final:  5,  Time:  1.8 
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Vortex  Filaments,  C, final:  10,  Time:  1.8 


3  4  5 

Vortex  Filaments,  C,final:  15,  Time:  1.8 


3.6  Model  Optimization 

The  simulation  is  optimized  by  varying  such  factors  as  the  number  of  vortex 
filaments,  the  spacing  between  vortex  filament  nodes,  the  time  step  in  the  simulation,  and 
the  location  of  the  vorticity  gradient  with  respect  to  the  numerical  end  of  the  vortex.  A 
balance  must  be  made  in  terms  of  resolution  of  the  model  and  numerical  efficiency. 

The  number  of  filaments  in  the  simulation  are  chosen  to  provide  an  initial 
vorticity  profile  that  very  nearly  matches  the  ideal  Gaussian  distribution,  while  at  the 
same  time  being  practical  to  implement  in  terms  of  computational  time.  Reference 
section  2.4. 1 .3  for  a  discussion  on  the  distribution  of  vorticity.  In  terms  of  balancing  the 
simulation  resolutions  with  simulation  runtime,  the  optimum  value  number  of  filaments  is 
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sixteen  (16)  vortex  filaments  evenly  distributed  about  the  circumference  of  nine  (9) 
radial  sleeves  of  vorticity.  Thus,  in  results  presented  here,  the  vortex  tube  is  discretized 
into  144  vortex  filaments. 

Next  the  length  of  the  vortex  filament  segment  is  optimized.  Again,  while  we 
would  like  the  segment  length  to  be  as  small  as  possible,  processing  constraints  require  a 
trade  off  between  small  segment  length  and  simulation  run  time.  Figure  3-1 8  shows  the 
axial  velocity  along  the  vortex  tube  centerline  at  T  =  1 .8  for  various  vortex  filament 
segment  lengths,  dz.  With  a  vortex  segment  length  of  dz  =  0.05,  the  simulation  run  time 
is  objectionable,  with  current  processing  power.  Therefore,  dz  =  0.1  is  used  in  the  results 
presented  in  the  previous  sections.  With  dz  =  0. 1,  the  simulation  exhibits  the 
characteristics  of  dz  =  0.05,  while  maintaining  a  practical  simulation  run-time.  For  dz  = 
0.2,  the  results  do  not  compare  well  with  the  dz  =  0.05,  so  this  vortex  segment  length  is 
not  used.  Therefore,  while  the  simulation  exhibits  qualitatative  convergence,  in  terms  of 
solution  characteristics,  it  does  not  exhibit  quantitative  convergence.  That  is  the  main 
reason  why  results  are  discussed  and  compared  in  terms  of  their  quality  and  not  the  exact 
numerical  quantity. 


BE 


0123456789 

Axial  dn-ection  (z) 

Figure  3-18  Simulation  optimization,  centerline,  axial  velocity,  vary  dz,  T  =  1.8 
We  can  see  this  qualitative  similarity  by  comparing  at  the  vortex  filament  plot  for 
these  cases.  Figure  3-19  shows  these  results  for  dz  =  0.2,  0.1,  and  0.05. 
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Figure  3-19  Vortex  filament  simulation  results,  vary  vortex  filament  segment  length,  dz 
The  same  optimization  process  is  done  with  the  simulation  time  step.  Figure  3-20 
shows  the  centerline  axial  velocity  with  the  time  step,  dt,  varying  from  0.1  to  0.001  at  T  = 
1.8.  Again,  to  achieve  a  balance  between  run-time  and  accuracy,  a  reasonably  small  time 
step  is  chosen.  This  reasonable  time  step  is  one  that  has  the  general  characteristics  of 
results  at  smaller  time  steps,  while  being  computational  efficient.  Thus,  for  this 


simulation,  we  choose  dt  =  0.01.  If  we  compare  the  vortex  filament  plots  that  results 
for  using  these  various  time  steps,  we  can  see  that  dt  =  0.01  is  a  reasonable  choice. 
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Axial  direction  (z) 

Figure  3-20  Simulation  optimization,  centerline  axial  velocity,  vary  dt,  T  =  1.8 
Figure  3-21  shows  the  vortex  filaments  at  T  =1.8  for  the  various  time  steps. 
Clearly,  when  dt  =  0.01,  the  results  take  on  the  characteristics  of  the  results  with  smaller 
time  step. 

Vortex  Filaments,  Time:  1 .8,  dt:  0.1 ,  dz:  0.1 
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Vortex  Filaments,  Time:  1 .8,  dt:  0.05,  dz:  0.1 
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Vortex  Filaments,  Time:  1 .8,  dt:  0.001 ,  dz:  0.1 


Figure  3-21  Vortex  filament  simulation  results,  vary  simulation  time  step,  dt 
Another  thing  that  is  optimized  is  the  location  of  the  end  of  the  vortex  tube. 
Theoretically,  the  vortex  tube  extends  downstream  to  infinity.  Numerically,  it  must  have 
some  finite  length.  However,  we  want  the  vortex  tube  to  end  far  enough  downstream  that 
that  end  has  a  negligible  effect  on  the  results.  The  area  of  interest  is  the  vorticity  gradient 
region,  thus  the  end  of  the  vortex  tube  is  measured  from  the  beginning  of  that  region. 
Figure  3-22  shows  the  velocity  along  the  centerline  of  the  vortex  tube  when  the  end  of  the 
vortex  tube  is  three  and  six  length  units  from  the  beginning  of  the  vorticity  gradient. 

There  is  virtually  no  change  in  the  result  when  the  length  to  the  numerical  end  of  the 
vortex  tube  is  doubled,  so  it  is  reasonable  to  use  a  length  of  three.  This  balances  the  need 
to  have  a  virtually  infinite  tube,  with  the  need  for  reasonable  computational  times.  Thus, 
in  the  results  previously  presented,  the  end  of  the  vortex  tube  is  three  length  umts  from 
the  beginning  of  the  vorticity  gradient  region.  In  temporal  terms,  with  a  freestream 
velocity  of  five  length  units/time  units,  the  vorticity  gradient  does  not  start  until  0.6  time 


units. 


no 


4.  Self-induction  theory  of  vortex  breakdown 

The  vortex  filament  method,  as  presented  in  the  previous  section,  indicates  vortex 
breakdown  behavior  in  the  vortex  tube.  The  only  mechanism  involved  in  that  model  is 
the  initial  vorticity  gradient  and  the  ensuing  self-induction  of  the  vortex  tube.  Thus,  we 
can  say  that,  if  the  vortex  filament  model  captures  vortex  breakdown,  the  cause  of  the 
vortex  breakdown  is  a  vorticity  gradient  and  the  mechanism  is  self-induction.  With  this 
self-induction  theory,  we  can  look  specifically  at  the  steps  that  lead  to  vortex  breakdown. 
Once  we  understand  these  steps,  we  can  better  understand  the  best  ways  to  prevent, 
delay,  or  weaken  vortex  breakdown. 

4.1  Difference  from  previous  theories 

Much  work  has  been  done  on  vortex  breakdown.  Many  of  these  ideas  focus  on 
the  shock-like  nature  of  vortex  breakdown,  with  many  people  making  the  analogy 
between  vortex  breakdown  and  a  hydraulic  jump.  Benjamin,  Keller,  Escudier,  Saffinan, 
and  others  have  down  much  exceptional  work  along  this  line  of  study.  We  do  not 
disagree  with  these  theories.  In  feet,  if  you  look  at  many  of  the  later  time  results  in  the 
vortex  filament  simulation,  you  can  see  many  shock-like  behaviors.  For  example,  you 
will  see  a  sudden  change  in  axial  velocity  down  the  centerline  of  the  vortex  tube.  Also, 
you  can  see  sudden  changes  in  the  sign  of  the  vorticity.  However,  our  model  does  not 
model  the  wave  effects  necessary  to  generate  shock  like  behavior.  Instead,  we  are 
modeling  the  unsteady  initial  formation  of  vortex  breakdown.  Thus,  once  a  steady  state 
breakdown  is  achieved,  modeling  the  vortex  breakdown  as  shock  or  hydraulic  jump  is 
likely  a  very  good  model.  However,  that  simplification  tells  you  little  about  the  actual 
formation  of  that  hydraulic  jump.  For  that,  we  find  that  the  inviscid  self-induction  theory 
provides  many  insights  not  previously  available.  Thus,  we  do  not  seek  to  supercede 
existing  vortex  breakdown  theory;  we  merely  bridge  the  gap  during  the  time  of  vortex 
breakdown  formation  and  the  steady  state  vortex  breakdown.  This  would  be  analogous  to 
studying  a  steady  shock  wave  with  simple  shock  relationships  versus  studying  the 
unsteady,  and  more  complicated  transient  formation  of  that  shock  wave. 
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4.2  Assumptions 

If  we  stai  achieve  vortex  breakdown  after  excluding  certain  physical  elements 
from  the  model,  we  can  i)e^ay  with  a  reasonable  amount  of  certainty  that  those  physical 
mechanisms  are  not  necessary  for  the  formation  of  vortex  breakdown.  In  this  case,  we 
are  limited  by  the  assumptions  implicit  in  the  vortex  filament  method.  Since  we  still 
achieve  vortex  breakdown  behavior,  we  can  say  that  the  physical  mechanisms  present  in 
the  vortex  filament  method  are  the  ones  important  to  the  vortex  breakdown  mechanism. 

In  the  vortex  filament  method,  we  have  assumed  unsteady,  incompressible,  inviscid  flow. 
The  motion  of  the  vortex  tube  is  determined  by  its  vorticity  distribution  and  the  induction 
of  the  vortex  tube  on  itself. 

4.3  Initial  cause  of  breakdown 

It  is  clear  from  the  results  that  the  initiator  of  vortex  breakdown  is  a  negative 
vorticity  gradient.  That  is,  the  vorticity  in  a  vortex  tube  is  greater  upstream  than 
downstream.  That  vorticity  gradient  will  be  in  the  axial  and  azimuthal  conponents  of 
vorticity,  though  any  axial  vorticity  gradient  will  cause  an  azimuthal  vorticity  gradient,  as 
Figure  4-1  shows.  In  Figure  4-1  we  see  that  the  axial  vorticity  will  induce  azimuthal 
velocities  in  the  vortex  tube.  These  azimuthal  velocities  with  twist  the  vortex  filaments. 
Because  the  vorticity  vector  is  always  aligned  with  these  vortex  filaments,  the  twisting  of 
the  vortex  filaments  results  in  a  component  of  vorticity  in  the  azimuthal  direction.  The 
effect  of  this  is  that  the  vortex  tube  develops  whatever  helical  angle  is  most  “natural”  to 
the  initial  distribution  of  axial  vorticity.  It  does  this  through  self-induction. 


Figure  4-1  Creating  Azimuthal  Vorticity 
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It  seems  clear  that  the  vorticity  gradient  is  the  prime  factor  in  initiating  the 
vortex  breakdown.  However,  we  can  take  it  a  step  further  and  ask,  what  causes  the 
vorticity  gradient?  If  the  cause  of  vortex  breakdown  is  self-induction  of  a  vortex  tube 
associated  with  a  negative  vorticity  gradient,  then  it  may  be  that  this  negative  vorticity 
gradient  can  be  caused  in  more  than  one  way.  Still,  the  negative  vorticity  gradient  is  the 
most  important  initiator  of  vortex  breakdown. 

4.3.1  Circulation  gradient 

In  the  vortex  filament  method  results  we  have  shown,  the  vorticity  gradient  is 
caused  by  a  circulation  gradient.  Vortex  filaments  are  linear  superimposed,  and  turn 
outward  at  the  appropriate  axial  location  to  preserve  solenoidality.  The  result  is  that  the 
upstream  region  of  vortex  filaments  has  more  vorticity  than  the  downstream  regions. 

This  creates  the  vorticity  gradient  that  leads  to  vortex  breakdown.  In  physical  situations, 
this  would  be  equivalent  to  a  pitch  up  and  hold  maneuver  over  a  delta  wing.  It  woiold 
also  be  equivalent  to  an  increase  in  velocity  over  a  delta  wing  at  a  fixed  angle  of  attack, 
since  the  increased  velocity  causes  increased  circulation.  In  a  vortex  generated  in  a  tube, 
this  would  be  equivalent  to  an  increase  in  pitch  of  the  swirl  vanes.  This  type  of  vorticity 
gradient  is  the  most  clear  and  straightforward  one  to  model  with  the  vortex  filament 
method.  That  is  why  we  limit  our  test  cases  to  this  type  vorticity  gradient  initiator. 

4.3.2  Velocity  gradient 

Over  strakes  and  delta  wing  aircraft,  one  can  think  of  another  way  to  achieve  this 
a  vorticity  gradient  without  an  explicitly  defined  circulation  gradient.  This  involves  the 
change  in  vorticity  associated  with  stretching  and  contraction  (the  opposite  of  stretching) 
in  a  vortex  tube  with  fixed  circulation.  Over  a  two-dimensional  wing,  there  is  a  peak 
negative  pressure  occurring  on  the  upper  surface  near  the  leading  edge  of  the  airfoil. 
From  this  peak,  the  pressure  increases  to  the  trailing  edge  of  the  airfoil,  with  positive 
pressure  occurring  at  the  trailing  edge.  Figure  4-2  shows  an  example  of  this  from  a 
classical  and  numerical  solution  for  the  pressure  coefficient  over  a  NACA  0012  airfoil  at 
9  degrees  angle  of  attack.  What  we  have  is  a  positive  pressure  gradient  (in  terms  of  the 


absolute  value  of  the  pressure  coefficient).  If  we  are  dealing  with  incompressible  or 
nearly  incompressible  flows,  we  can  talk  equivalently  about  velocity,  since  the  pressure 
coefficient  is  defined  by  equation  (1-55). 


So  if  we  have  a  positive  pressure  gradient,  we  also  have  a  negative  velocity 
gradient.  Thus,  if  you  were  to  place  a  vortex  tube  in  the  pressure  field,  it  would  see  a 
velocity  that  decreases  along  the  downstream  axial  direction  of  the  vortex  tube.  Figure 
4-2  shows  the  velocity  that  the  vortex  tube  sees  downstream  of  the  negative  pressure 
peak. 


Figure  4-2  Pressure  distribution  over  NACA  0012  airfoil  {ref.  5) 
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Figure  4-3  Vortex  tube  in  negative  velocity  gradient 
This  negative  velocity  gradient  will  affect  the  vorticity  in  a  vortex  tube  in  terms  of 
stretching.  The  vorticity  equation  defines  this  interaction,  as  in  equation  (4-2).  We  will 
only  take  the  stretching  term  into  account  and  if  we  will  assume  velocity  and  vorticity 
only  change  in  the  axial  direction  (z-direction).  Also,  we  will  assume  that  the  velocity  is 
constant,  while  the  vorticity  can  change  in  time.  If  we  do  this  we  get  the  expression  in 
equation  (4-3). 


Dt  ^ 

(4-2) 

dt  dz  dz 

(4-3) 

If  we  were  looking  at  a  steady  state  case,  the  time  changing  term  would  drop  out, 
and  we  would  have  equation  (4-4). 

(«) 

^  ^  dz  "  & 

If  we  assume  that  u{z)  and  &^{t)  remain  the  same  sign  (aligned  in  the  axial  direction), 

we  can  very  easily  see  that  the  and  terms  will  have  to  be  the  same  sign.  If 

dz  dz 

the  positive  z-direction  is  to  the  right  (downstream),  then  a  negative  velocity  gradient  will 
cause  the  vortex  tube  to  have  a  negative  vorticity  gradient.  This  is  similar  to  the 
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argument  that  M. J.  Lighthill  introduces  and  Whitman  derives  when  they  show  that 
vorticity  vectors  not  only  move  with  the  fluid  elements,  but  are  proportional  to  the  length 

of  the  line  element,  8/  {ref.  23).  Therefore,  c5  =  c  •  8f ,  so  if  velocity  changes  along  the 
wing,  the  length  of  line  element,  8/ ,  will  change,  as  will  the  magnitude  of  the  vorticity 
vector.  Thus,  it  is  clear  that  the  velocity  profile  above  wing  associated  with  an  adverse 
pressure  gradient  is  sufficient  to  cause  a  vorticity  gradient  in  a  vortex  tube.  With  that 
vorticity  gradient,  this  vortex  tube  can  cause  itself  to  breakdown  through  the  self- 
induction  mechanism. 

4.4  Self-induced feedback  mechanism 

If  a  vortex  tube  causes  itself  to  breakdown,  what  steps  does  it  take  to  do  so? 
Others  have  proposed  various  feedback  mechanisms  for  vortex  breakdown  {ref  20),  but 
we  will  approach  it  fi-om  a  self-induction  perspective,  with  the  starting  point  being  the 
vorticity  gradient,  however  it  is  brought  about.  With  this,  we  have  the  following 
feedback  mechanism  that  causes  vortex  breakdown.  Using  the  results  Ifrom  the  vortex 
filament  simulation,  we  will  explain  each  part  of  this  diagram  in  this  section. 
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Figure  4-4  Self-Induced,  inviscid  feedback  mechanism  of  vortex  breakdown 
4.4.1  Radial  expansion 

As  shown  in  previous  sections,  an  axial  vorticity  gradient  will  lead  to  an 
azimuthal  vorticity  gradient.  If  we  look  at  a  cross  section  of  a  vortex  tube  in  the  region  of 
an  aTimiitbal  vorticity  gradient,  we  will  see  what  is  illustrated  in  Figure  4-5.  The  red 
circles  indicate  azimuthal  vorticity,  ©e,  out  of  the  plane.  The  blue  crosses  indicate 
azimuthal  vorticity  into  the  plane.  Their  arced  arrows  show  the  individual  induced 
velocities.  The  black  arrows  show  the  net  induced  radial  velocity  along  the  vortex  tube. 
The  azimuthal  vorticity  increases  in  the  upstream  direction  to  the  left.  In  this  situation, 
we  can  see  that  the  induced  radial  velocity  is  greater  downstream  than  upstream. 
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Figure  4-5  Radial  expansion 

This  becomes  clearer  if  we  look  at  the  induced  radial  velocities  in  more  detail. 
Figure  4-6  shows  an  expanded  diagram  of  a  cross  section  of  the  vortex  tube  in  the  region 
of  the  vorticity  gradient.  We  will  represent  the  azimuthal  vorticity  component  as  2-D 
potential  vortices,  with  each  vortex  spaced  a  distance,  d,  apart.  The  strength  of  the 
vorticity  is  Tj.  In  this  simplified,  2-D  case,  the  induced  velocity  from  one  vortex  on 


Figure  4-6  Azimuthal  vorticity  in  meridional  plane 
The  vorticity  gradient  can  be  represented  by  increasing  T  to  the  left.  Table  4-1 
shows  the  induced  radial  velocities  associated  with  this  setup. 
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Table  4-1  Induced  radial  velocities 
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Figure  4-7  plots  the  induced  radial  velocities  from  Table  4-1.  Looking  at  even 
this  simple  case,  we  can  see  how  a  vorticity  gradient  induces  radial  velocities  that  will 
lead  to  radial  expansion.  The  large  negative  induced  velocity  for  vortex  8  is  largely 
caused  by  vortex  7  and  the  absence  of  a  vortex  to  the  left  of  vortex  8.  In  essence,  this  is 
similar  to  the  radial  expansion  and  contraction  between  two  vortex  rings,  as  in  the  case  of 
the  so-called  “leap-frogging”  vortex  rings. 
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Figure  4-7  Induced  radial  velocities 

Now  we  can  expand  this  2-D  analogy  to  the  3-D  case.  As  a  result  of  this  induced 
radial  velocity  caused  by  the  azimuthal  vorticity  gradient,  the  vortex  tube  expands 
radially  in  the  vorticity  gradient  region.  This  is  exactly  what  happens  in  the  3-D  vortex 
filament  simulation.  Figure  4-8  shows  the  rendered  vortex  tube  for  the  filaments  initially 
at  r  =  0.12  and  r  =  0.42.  In  both  cases,  by  T  =  0.9  (the  increase  in  vorticity  started  at  T  = 
0.6),  the  vortex  tube  has  started  to  expand  radially.  The  red  shaded  region  shows  the  area 
where  the  circulation  is  increasing.  The  radial  spokes  of  vorticity  that  maintain 
solenoidaUty  of  the  vorticity  field  are  left  out  of  these  images  for  clarity.  The  radial 
expansion  continues.  However,  it  also  leads  to  several  other  effects,  like  pile-up  and 
vorticity  sign  switching  that  are  visible  in  this  image.  We  cannot  isolate  these  effects  of 
the  feedback  mechanism,  because  they  all  interact  with  each  other  simultaneously  to  lead 
to  vortex  breakdown. 
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Figure  4-8  Radial  expansion,  rendered  vortex  filaments,  r  =  0.12  (L),  r  =  0.42  (R),  T  =  0.9 
If  we  look  at  the  radial  velocity  component  at  this  time,  we  should  see  the 
velocities  that  lead  to  this.  Figure  4-9  shows  the  color  contour  plot  of  the  induced  radial 
velocity  for  this  area.  The  induced  outward  radial  velocity  is  clearly  visible  in  the 
vorticity  gradient  region. 
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4.4.2  Pile-up 

This  radial  expansion  leads  to  the  situation  we  call  pile-up.  Pile-up  refers  to  the 
axial  contraction  of  the  vortex  tube.  However,  it  should  be  noted  that,  even  if  the  radial 
expansion  were  suppressed  somehow,  pile-up  would  still  occur  in  regions  of  an  azunuthal 
vorticity  gradient.  We  will  first  ignore  the  radial  expansion,  and  then  we  will  take  it  into 
account,  showing  how  it  increases  the  effect  of  pile-up.  Pile-up  happens  as  a  result  of  an 
azimuthal  vorticity  gradient.  The  azimuthal  vorticity  gradient  causes  an  axial  velocity 
gradient,  as  shown  in  Figure  4-10.  The  red  circles  indicate  azimuthal  vorticity,  coe,  out  of 
the  plane.  The  blue  crosses  indicate  azimuthal  vorticity  into  the  plane.  Their  arced 
arrows  show  the  individual  induced  velocities.  The  black  arrows  show  the  net  induced 
axial  velocity  along  the  centerline  of  the  vortex  tube.  The  azimuthal  vorticity  increases  in 
the  upstream  direction  to  the  left.  Like  a  pile-up  with  automobiles  on  a  highway,  the 
vortex  tube  upstream  is  moving  slower  in  the  axial  direction  than  the  vortex  tube 
upstream. 
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Figure  4-10  Vorticity  gradient  leading  to  pile-up 
This  axial  velocity  gradient  causes  the  axial  components  of  the  vortex  filaments  to 
pile-up.  When  the  axial  vorticity  components  pile  up,  they  undergo  contraction  (the 
opposite  of  stretching).  This  contraction  leads  to  a  further  decrease  is  axial  vorticity  in 
the  pile-up  region.  Thus,  a  azimuthal  vorticity  gradient  in  the  axial  direction  causes  pile- 
up  and  a  decrease  in  axial  vorticity  in  that  region.  The  decrease  in  axial  vorticity  leads  to 
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a  corresponding  decrease  in  azimuthal  vorticity,  which  leads  to  a  stronger  azimuthal 
vorticity  gradient  and  more  pile-up.  Due  to  this  feedback  response,  the  pile-up  effect  is 
magnified  over  time. 

Not  only  are  the  two  effects  of  radial  e}q)ansion  and  pile-up  coupled,  they  are  also 
accentuated  by  each  other.  Figure  4-1 1  shows  a  diagram  of  this. 
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Figure  4-11  Pile-up  with  radial  expansion 
The  reason  a  radially  expanded  vortex  tube  increases  the  effect  of  pile-up  is  that 
the  a^imiithfll  voiticity  con^onents  have  expanded  radially  and  moved  farther  away  from 
the  vortex  tube  centerline.  As  a  result,  these  components  induce  less  axial  velocity  along 
the  centerline.  As  a  result  of  radial  expansion,  induced  velocity  along  the  centerline  of 
the  vortex  tube  slows  down  even  more.  The  upstream  axial  velocities  do  not  change,  so 
the  effects  of  pile-up  increase.  This  is  one  way  in  which  the  feedback  mechanisms  are 
linked. 

We  can  see  this  effect  of  the  vorticity  gradient  and  the  radial  expansion  in  the 
vortex  filament  simulation.  The  first  figure,  Figure  4-12,  shows  the  vortex  tube  at  T  = 

0.9  and  T  =  1.2.  Here,  the  red  line  shows  tracks  particles  along  the  centerline  of  the 
vortex  tube.  The  cross  marks  show  distance  between  the  tracking  line  segment  endpoint 
nodes.  They  were  originally  spaced  at  0.1  length  umts  apart.  The  effect  of  pile-up  is 
clear  visible  in  the  z  =  1 .5  to  2.0  region,  which  is  the  region  of  the  vorticity  gradient. 
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Vortex  Filaments,  Time:  0.9 


Figure  4-12  Pile-up  along  vortex  tube  centerline  (top)  T  =  0.9,  (bottom)  T  =  1.2 
We  can  see  how  this  fits  in  to  the  overall  vortex  tube  by  overlaying  the  r  =  0.12 
filaments  on  the  centerline  plot  in  Figure  4-13.  Clearly,  the  pile-up  is  strongest  in  the 
region  of  the  largest  radial  expansion. 


Vortex  Filaments,  Time:  0.9 


Vortex  Filaments,  Time:  1 .2 


Figure  4-13  Pile-up  along  vortex  tube,  center  and  r  =  0.12,  (top)  T  ~  0.9,  (bottom)  T  -  1.2 
If  we  look  at  the  rendered  view  of  the  same  radius,  in  Figure  4-14,  we  can  see  that 
pile-up  occurs  throughout  the  vortex  tube  and  not  just  at  the  centerline.  Looking  at  the 
circumferential  tracking  lines,  we  see  that  the  vortex  filament  segments  endpoint  nodes 
get  closer  and  closer  together  and  they  get  nearer  to  the  vorticity  gradient  region.  This  is 
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a  result  of  pile-up.  As  a  consequence,  we  should  expect  the  axial  component  of 
vorticity  to  decreases  throughout  the  vortex  tube  in  regions  of  pile-up. 

Vortex  Filaments,  Time:  0.9 
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Figure  4-14  Pile-up  along  vortex  tube,  rendered  r  =  0.12,  (top)  T  =  0.9,  (bottom)  T  =  1.2 
The  effects  of  pile-up  will  be  clear  in  the  plots  of  axial  velocity  along  the 
meridional  plane.  We  should  also  be  able  to  see  the  decrease  in  axial  vorticity  on  this 
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Axial  Vorticity,  meridional  plane,  Time:  0.9 
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Axial  Vorticity,  meridional  plane.  Time:  1 .2 


Figure  4-16  Pile-up,  axial  vorticity,  (top)  T  =  0.9,  (bottom)  T  =  1.2 
Though  we  have  tried  to  talk  about  radial  expansion  and  pile-up  separately,  they 
are  really  inseparable.  One  inevitably  leads  to  the  other  due  to  continuity.  The  continuity 
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equation  for  the  vortex  tube  is  — —  J+ 

r  or  oz 


=  0 .  Due  to  the  negative  azimuthal 


du. 


vorticity  gradient,  there  is  a  negative  velocity  gradient,  and  pile-up,  so 


—iru  )  >  0  in  order  for  the  continuity  equation  to  be  satisfied.  If  we  integrate  —  {ru^ ) 
dr  ^ 


from  0  to  fl  where  a  is  the  radius  of  the  vortex  tube,  we  find  that 
a-uXr  =  a)-auXr  =  0)>^-  Of  course,  M^(r  =  0)=0,so  a-uXr  =  a)>Q  must  be 
satisfied.  The  expression  «,.  >  0 ,  implies  radial  eiqiansion.  Thus,  through  continuity, 
radial  expansion  and  pile-up  are  inseparable. 


4. 4. 3  Vorticity  sign  switch 

As  a  result  of  pile-up  and  radial  expansion,  there  is  a  decrease  in  axial  vorticity  in 
the  region  of  the  vorticity  gradient.  Axial  vorticity  is  associated  with  azimuthal  velocity. 
Therefore,  in  the  regions  of  pile-up,  radial  expansion,  and  the  axial  vorticity  gradient,  the 
vortex  tube  will  rotate  more  slowly  than  the  upstream  region  of  the  vortex  tube.  You  can 
also  think  of  this  in  terms  of  conservation  of  angular  momentum.  If  a  vortex  tube 
increases  in  radius,  it  will  have  to  rotate  more  slowly  in  that  region  to  conserver  angular 
momentum.  This  decrease  in  azimuthal  velocity  in  the  region  of  radial  e3q)ansion  and 
pile-up  is  clearly  visible  in  the  color  contour  plot  of  this  velocity  along  the  meridional 
plane.  It  is  shown  in  Figure  4-17. 


Azimuthal  Velocity,  meridional  plane,  Time:  0.9 
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This  azimuthal  velocity  gradient  causes  tilting  in  the  vortex  filaments  in  the 
region  of  radial  expansion  and  pile-up.  The  vorticity  vector  is  at  all  times  aligned  with 
the  vortex  filament.  Therefore,  due  to  this  tilting,  the  vector  decomposition  of  the  vortex 
filament  will  show  a  sign  switch  in  the  azimuthal  vorticity.  Brown  and  Lopez  were  the 
first  to  mention  this  very  important  result  {ref.  1 7).  We  confirm  that  result  here.  To 
explain  this  effect,  let  us  look  at  the  diagram  in  Figure  4-18. 


Figure  4-18  Vortex  tube,  rotating  and  different  velocities 
Here,  we  show  the  effect  of  an  azimuthal  velocity  gradient  on  a  portion  of  a 
vortex  tube.  The  vortex  tube  is  fixed  at  its  source.  At  a  location,  P,  downstream  of  that 
source,  it  is  rotating  at  a  speed  ue,i.  Downstream  of  P,  the  vortex  tube  has  undergone 
pile-up  and  radial  expansion.  Thus  at  Q,  it  rotates  at  a  slower  speed,  ue,2.  Now,  let  us 
track  a  vortex  filament  in  this  vortex  tube.  For  simplicity,  let  us  say  that  at  one  instant,  to, 
the  vortex  filament  is  axially  aligned  with  the  vortex  tube.  We  will  track  this  the  motion 
of  this  vortex  filament  at  subsequent  times.  Figure  4-19  shows  how  it  will  move.  In  this 
figure,  the  vortex  line  is  green.  Because  it  is  fixed  at  the  source  and  because  it  rotates 
faster,  as  indicated  by  ue.i,  we  can  see  that  the  vortex  filament,  which  moves  with  the 
fluid,  will  develop  a  “v-shaped”  kink. 
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Figure  4-19  Azimuthal  vorticity  sign  switch 
The  vorticity  vector  is  always  aligned  with  this  vortex  filament,  so  if  we 
decompose  the  vector  components  of  the  vortex  filament,  we  can  see  how  the  vorticity 
changes.  The  axial  component  of  vorticity  does  not  change,  but,  as  the  above  figure 
shows,  the  azimuthal  vorticity,  a)e»  undergoes  a  sign  switch.  This  is  precisely  what 
happens  in  the  vortex  filament  simulation.  We  can  see  this  by  looking  at  a  single  vortex 
filament  fi'om  the  vortex  filament  simulation  at  T  =  1 .2.  Figure  4-20  shows  the  vortex 
filaments  that  were  originally  at  r  =0.12. 
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Figure  4-20  Vortex  filaments,  azimuthal  vorticity  sign  switch 
In  this  view  we  see  the  vortex  filaments.  The  yellow  line  highlights  one  specific 
filament.  If  we  decompose  that  filament  into  its  components,  we  can  see,  as  was 
illustrated  in  Figure  4-19,  the  azimuthal  vorticity  sign  switch.  This  vorticity  sign  switch 
is  especially  clear  in  the  vorticity  that  is  calculated  along  the  meridional  plane.  Figure 
4-21  shows  the  azimuthal  vorticity  in  the  meridional  plane  at  T  =  1.2. 
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Figure  4-21  Azimuthal  vorticity  sign  switch,  meridional  plane 
This  azimuthal  vorticity  sign  switch  is  important  to  vortex  breakdown  because  it 
is  the  condition  that  causes  backflow  within  this  vortex  tube.  Without  backflow,  vortex 
breakdown  would  be  a  transient  phenomenon  instead  of  one  that  can  exist  at  a  fixed 
location  relative  to  an  aircraft  wing  in  steady  state. 

We  can  understand  how  the  azimuthal  vorticity  sign  switch  leads  to  backflow  by 
looking  at  the  diagram  in  Figure  4-22.  Instead  of  inducing  axial  velocities  in  the 
direction  of  the  fi-eestream,  the  region  where  the  azimuthal  vorticity  has  switched  its  sign 
induces  velocities  pointing  in  the  upstream  direction.  The  backflow  becomes  stronger  as 
the  ayimiitbal  vorticity  sign  switching  becomes  stronger.  When  the  induced  backflow 
equals  the  relative  freestream  velocity,  the  vortex  breakdown  region  will  stop  moving 


downstream. 
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Figure  4-22  Backflow 

We  can  see  the  start  of  this  backflow  in  the  axial  velocity  plot  along  the  meridional  plane, 
as  in  Figure  4-23.  We  previously  looked  at  this  figure  in  terms  of  pile-up. 


Axial  Velocity,  merldonal  plane,  Time:  1 ,2 


Figure  4-23  Beginning  of  backflow,  axial  velocity,  T  =  1.2 
We  can  now  see  how  these  mechanisms  enhance  one  another  in  a  feedback  loop. 
As  the  pile-up  continues,  the  vortex  radially  expands  and  the  azimuthal  vorticity  switches 
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signs  to  induce  backflow.  This  backflow  enhances  the  pile-up,  which  enhances  the 
radial  expansion,  which  enhances  the  vorticity  sign  switch,  which  in  turn  enhances  the 
backflow. 

It  can  also  be  noted  that  a  radial  vorticity  sign  switch  is  associated  with  radial 
expansion.  This  is  because  the  vorticity  vector  is  aligned  at  all  times  with  the  vortex 
filaments.  Therefore,  in  regions  of  radial  expansion,  there  will  be  a  positive  radial 
component  of  vorticity  just  as  the  radius  of  the  vortex  tube  increases.  Likewise,  there 
will  be  a  negative  radial  component  of  where  the  vortex  tube  contracts.  The  position  of 
zero  radial  vorticity  (co,.  =  0  )  corresponds  to  the  maximum  radii  of  the  vortex  filaments. 
Here,  the  radial  expansion  ends  and  radial  contraction  begins.  Beyond  this  point,  the 
vortex  tube  starts  to  close  on  the  downstream  end,  as  compared  to  the  hypothesized  open 
cavity  model  adopted  by  Keller  et  al  (see  ref.  24  and  25  for  details  on  the  open  cavity 
model).  Figure  4-24  shows  the  resulting  radial  vorticity  sign  switch. 


Radial  Vorticity,  meridonal  plane,  Time:  1 .2 


Figure  4-24  Radial  vorticity  sign  switch 
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4.4.4  Turning  point 

The  azimuthal  vorticity  sign  switch  induces  velocities  that  cause  backflow  inside 
the  vortex  tube.  However,  this  backflow  also  causes  vortex  filaments  to  turn  inward  on 
themselves.  We  can  see  how  this  happens  by  looking  at  the  illustration  in  Figure  4-25. 
The  azimuthal  vorticity  sign  switch  induces  a  turning  of  the  vortex  filaments.  The 
filaments  turn  inward,  leading  to  a  radial  contraction  of  the  vortex  tube  as  in  Figure  4-25 
(L).  However,  they  do  not  stop  turning  at  this  point.  They  continue  to  turn,  being 
brought  upstream,  relative  to  the  region  of  the  azimuthal  vorticity  sign  switch,  by  the 
induced  backflow  velocities  as  in  Figure  4-25  (R).  This  radial  contraction  is  important 
because,  when  combined  with  the  backflow,  it  leads  to  the  culminating  point  of  the 
inviscid  feedback  mechanism— the  turning  point. 


Figure  4-25  Turning  point  formation 

There  is  a  definite  location  where  the  vortex  filaments  have  turned  inward  on 
themselves,  as  in  Figure  4-26.  We  call  this  the  turning  point.  The  exact  location  of  the 
turning  point  is  defined  as  where  the  vortex  filament  is  vertical  relative  to  the  axial 
direction. 
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Figure  4-26  Turning  point  established 

At  this  point,  the  axial  vorticity  is  zero.  Upstream  of  the  turning  point,  along  the 
vortex  filament,  the  axial  vorticity  is  positive.  Downstream  of  the  turning  point,  along 
the  vortex  filament,  the  axial  vorticity  is  negative.  It  becomes  positive  again  near  the 
centerline  of  the  vortex  tube  as  the  vortex  filaments  turn  back  downstream. 

We  can  also  see  this  turning  point  most  clearly  in  the  vortex  filament  results.  If 
we  look  at  a  projection  of  the  vortex  filaments  onto  the  meridional  half-plane,  we  can 
clearly  see  the  location  of  the  turning  point  at  various  radii  of  the  vortex  tube.  Figure 
4-27  shows  these  results.  The  location  of  the  turning  point  is  marked  with  a  vertical  Ime. 


Figure  4-27  Vortex  filaments  projected  on  meridional  plane,  T  =  1.5 
The  projection  of  the  vortex  filaments  onto  the  meridional  half-plane  allows  us  to 
see  the  direction  of  axial  and  radial  vorticity  vectors,  as  they  are  aligned  with  the  vortex 
filaments.  Thus,  we  can  see  the  axial  vorticity  sign  switch.  Figure  4-28  shows  the  axial 
and  radial  components  of  the  vortex  filaments  overlaid  with  the  rendered  vortex  filaments 
fi-om  the  r  =  0.12  initial  radius  at  T  =  1 .5. 


Vortex  Filaments,  Time:  1 .5 


A  >  tiirniiia  ooint 


Figure  4-28  Axial  vorticity  sign  switch 

We  can  also  see  the  axial  vorticity  sign  switch  in  the  data  along  the  meridional 
plane.  Figure  4-29  shows  this  at  T  =  1.5.  The  axial  vorticity  sign  switch  corresponds  to 
the  area  near  the  turning  point  throughout  the  vortex  tube. 
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Axial  Vorticity,  meridional  plane.  Time:  1 .5 


Figure  4-29  Axial  vorticity  sign  switch,  meridional  plane,  T  =  1.5 

This  axial  vorticity  sign  switch,  interior  to  the  vortex  tube,  induces  counter 
rotating  flows  interior  to  the  vortex  tube.  It  is  important  to  note  that  these  counter 
rotating  flows  exist  only  interior  to  the  vortex  tube  and  in  the  region  of  the  turning  point. 
They  do  not  exist  outside  of  the  vortex  tube.  These  counter  rotating  flows  lead  to  regions 
where  viscosity  becomes  important  with  shear  layers  and  the  development  of  small-scale 
turbulence.  As  a  result,  the  inviscid  model  reaches  the  limits  of  its  assumptions  after  the 
turning  point  forms.  However,  with  the  turning  point,  all  features  of  vortex  breakdown, 
as  compared  to  experimental  data,  are  present.  We  will  compare  the  numerical  results  to 
experimental  data  in  later  sections. 

In  summary,  vortex  breakdown  has  been  initiated  to  a  recognizable  form  through 
an  entirely  inviscid  mechanism  of  pile-up,  radial  expansion  and  contraction,  and  vorticity 
sign  switching.  The  vorticity  sign  switch  is  present  in  all  forms  of  vorticity  ^azimuthal, 
radial,  and  axial.  These  inviscid,  self-induced  mechanisms  work  together  to  form  the 
turning  point.  Once  the  turning  point  forms,  the  downstream  region  of  the  vortex 
breakdown  forms  a  dissipative,  turbulent  flow  that  leads  to  a  steady  state  vortex 
breakdown.  Thus,  vortex  breakdown  is  initiated  by  a  vorticity  gradient  that  manifests 
itself  through  and  inviscid,  self-induced  feedback  mechanism.  This  inviscid  feedback 
mechanism  continues  xmtil  it  forms  the  turning  point.  This  turning  point  and  its 
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associated  axial  vorticity  sign  switch  cause  counter  rotating  flows  interior  to  the 
vortex  tube,  which  leads  to  viscous  actions  and  steady  state  vortex  breakdown.  This 
leads  to  the  conclusion  that  if  we  can  suppress  a  portion  of  the  inviscid  feedback 
mechanism,  it  would  seem  possible  to  prevent,  weaken,  or  delay  the  fijrmation  of  the 
turning  point,  and  by  doing  so,  prevent,  weaken,  or  delay  the  formation  of  vortex 
breakdown. 

4.5  Inqtlications  of  self-induction  mechanism  on  vortex  breakdown 

If  we  understand  that  vortex  breakdown  works  through  a  self-induced,  inviscid 
mechanism,  we  can  better  understand  why  vortex  breakdown  occurs  at  specific 
chordwise  locations  over  a  delta  wing.  If  the  induced  backflow  velocity  is  caused  by 
azimuthal  vorticity  tilting,  it  should  follow  that  for  a  given  circulation  in  a  vortex  tube, 
there  is  a  mfivimimn  possible  induced  backflow,  depending  on  the  degree  of  the  azimuthal 
vorticity  sign  switch.  This  maximum  induced  backflow  is  a  function  of  the  circulation  of 
the  vortex  tube,  which  itself  is  a  function  (in  terms  of  a  delta  wing)  of  the  delta  wing 
sweep  angle,  the  freestream  velocity,  the  delta  wing  angle  of  attack,  and  yaw  angle.  The 
other  important  factor  is  the  velocity  profile  along  the  vortex  tube  path.  As  we  saw 
before,  that  velocity,  along  the  path  of  the  vortex  tube,  reaches  a  peak  near  an  airfoil 
leading  edge,  and  then  decreases  toward  the  trailing  edge.  This  velocity  profile,  itself,  is 
a  function  of  the  delta  wing  sweep  angle,  the  fi-eestream  velocity,  and  the  delta  wing 
angle  of  attack. 

Therefore,  for  a  given  delta  wing  sweep  angle,  fi-eestream  velocity,  delta  wing 
angle  of  attack,  and  yaw  angle,  there  is  only  one  location  over  the  wing  where  vortex 
breakdown  can  remain  fixed  with  reference  to  the  wing.  That  point  is  where  the  relative 
axial  fi-eestream  velocity  exactly  equals  the  maximum  induced  backflow  velocity  as  it 
manifests  itself  in  the  vortex  tube.  The  relative  axial  fireestream  velocity  is  the  fi-ee 
stream  velocity  that  the  region  of  vortex  breakdown  experiences.  This  is  relative  to  the 
position  of  the  vortex  breakdown  over  the  wing  since  the  fi-eestream  velocity  changes 
near  the  vortex  along  the  path  of  the  vortex  over  the  wing.  Figure  4-30  illustrates  this 
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idea,  where  point  P  corresponds  to  such  a  point  of  balance  between  the  maximum 
induced  backflow  and  the  relative  axial  freestream  velocity. 

I  relative  axial  freestream  velocity 


relative  axial  freestream  velocity 


Figure  4-30  Steady  state  vortex  breakdown  location 
This  steady  vortex  breakdown  location  would  represent  a  very  delicate  point, 
where  a  slight  change  in  induced  backflow  velocity  or  of  the  relative  axial  freestream 
velocity  would  require  the  vortex  breakdown  location  to  rebalance. 

With  this  idea  we  can  explain  how  the  vortex  breakdown  location  changes  as 
results  of  an  increase  in  angle  of  attack.  An  increase  in  angle  of  attack  will  increase  the 
circulation  in  the  vortex  tube  and  change  the  relative  axial  freestream  velocity  profile  by 
increasing  the  peak  velocity  and  decreasing  the  velocities  downstream  of  the  peak.  As  a 
result  of  the  increase  in  circulation,  the  maximum  induced  backflow  velocity  will 
increase.  Thus,  as  we  see  in  Figure  4-31,  the  steady  state  vortex  breakdown  will  move 
forward  from  P  to  P’. 
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Figure  4-3 1  Vortex  breakdown  location  change  with  AOA  increase 
A  decrease  in  angle  of  attack  does  just  the  opposite.  It  will  cause  the  vortex 
breakdown  location  to  move  downstream.  If  the  angle  of  attack  is  sufficiently  small,  the 
maximum  induced  backflow  will  not  be  sufficient  to  maintain  the  vortex  breakdown  over 
the  wing.  It  will  be  swept  downstream  where  it  will  dissipate.  The  result,  at  lower  angles 
of  attack,  is  no  vortex  breakdown  over  the  wing. 

A  change  in  freestream  velocity,  with  the  angle  of  attack  fixed,  results  in  little 
change  in  angle  of  attack  (provided  that  the  Reynolds  number  is  sufficiently  low)  because 
the  circulation  and  corresponding  maximum  induced  backflow  increases  proportionally 
to  the  amount  that  the  change  in  the  relative  axial  freestream  velocity  profile,  as  in  Figure 


4-32. 
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chordwise  location 


Figure  4-32  Vortex  breakdown  location  change  with  freestream  velocity  increase 
Typically,  the  maximum  induced  backflow  and  the  relative  axial  freestream 
velocity  are  coupled  through  the  freestream  velocity,  the  delta  wing  angle  of  attack,  yaw 
angle,  and  the  delta  wing  sweep  angle.  However,  by  understanding  how  they  work 
together  through  the  self-induction  mechanism,  we  can  come  up  with  various  schemes  to 
decouple  these  to  effects,  and  thus  provide  a  way  to  control  the  location  of  vortex 
breakdown.  For  example,  if  we  could  reduce  the  circulation  in  a  vortex  tube,  while 
keeping  the  relative  axial  freestream  velocity  fixed,  we  would  be  able  to  move  the  vortex 
breakdown  location  downstream.  This  type  of  vortex  breakdown  control  mechanism 
would  provide  a  way  to  use  the  location  of  the  vortex  breakdown  for  aircraft  control 


purposes. 
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5.  Experimental  and  Computational  Evidence 

While  the  vortex  filament  method  simulation  is  a  useful  tool  for  understanding 
vortex  breakdown  and  the  self-induction  theory  of  vortex  breakdown,  it  needs  to  be 
validated  with  experimental  data.  We  will  vaUdate  the  vortex  filament  simulation  results 
qualitatively  with  several  experiments.  The  first  experiment  was  performed  by  S. 
Srigrarom  specifically  to  match  the  vorticity  gradient  case.  See  ref.  4  for  conqilete 
details.  It  involves  a  delta  wing  at  fixed  angle  of  attack,  undergoing  a  transient  increase 
in  fi-eestream  velocity.  This  change  in  freestream  velocity  creates  a  vorticity  gradient  that 
leads  to  vortex  breakdown  {ref  4).  The  second  experiment  Ulustrates  how  vorticity 
gradients  in  vortex  tubes  of  fixed  circulation  can  lead  to  vortex  breakdown.  The  third 
case  was  a  previously  published  case  of  vortex  breakdown  in  a  cylinder  with  a  rotating 
end  wall  by  M.P.  Escudier  {ref  27).  This  case  was  chosen  because  of  its  similarity  to  the 
conditions  in  the  vortex  filament  simulation.  The  fourth  case  is  a  full  Navier-Stokes 
simulation  of  the  Escudier  eiqieriment  {ref.  28).  This  allows  us  to  compare  qualitatively 
the  velocity  and  vorticity  data  of  a  fully  viscous  solution  of  vortex  breakdown  with  the 
inviscid  solution  from  the  vortex  filament  simulation.  Finally,  previously  published 
experiments  involving  vortex  breakdown  control  will  be  examined  in  terms  of  the  self- 
induction  theory  of  vortex  breakdown. 

5.1  Case  I,  vortex  breakdown  over  delta  wing  with  vorticity  gradient 

To  validate  the  vortex  filament  simulation,  we  will  compare  the  vortex  filament 
simulation  to  the  vortex  breakdown  over  a  delta  wing  caused  by  a  vorticity  gradient  in  a 
water  tunnel.  This  experiment  is  different  than  the  simulation,  because  it  involves  a  wing 
rather  than  just  a  vortex  tube.  This  introduces  asymmetries  that  are  not  present  in  the 
simulation.  Because  of  these  differences  between  the  experiment  and  the  simulation,  we 
will  make  qualitative  and  general  comparisons,  rather  than  qiiantitative  conqiarisons.  In 
the  experiment  dye  flow  visualization,  laser  induced  fluorescence  (LIF),  and  particle 
im^e  velocimetry  (PIV)  are  used  to  analyze  the  flow. 
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5.1.1  Experiment  setup 

This  experiment  used  a  Plexiglas  delta  wing  with  a  60  degrees  apex  angle 
mounted  at  20  degrees  angle  of  attack  in  a  closed  circuit  water  tunnel.  A  vorticity 
gradient  was  generated  by  quickly  increasing  the  freestream  velocity  in  the  water  tunnel 
from  0  to  5  cm/s.  This  represents  a  Reynolds  number  change  from  0  to  20,000  based  on 
the  delta  wing  chord.  Figure  5-1  shows  how  this  change  in  freestream  velocity  creates  a 
negative  vorticity  gradient. 


Figure  5-1  Generating  a  vorticity  gradient  over  a  delta  wing 
This  negative  vorticity  gradient  is  what  leads  to  vortex  breakdown  in  the  vortex 
filament  simulation.  If  the  initial  cause  of  vortex  breakdown  is  a  negative  vorticity 
gradient,  then  the  vorticity  gradient  caused  by  the  freestream  velocity  increase  should 
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also  lead  to  vortex  breakdown  in  the  region  of  the  vorticity  gradient.  This  is  precisely 
what  is  seen  in  this  delta  wing  experiment. 

5.1.2  Dye  flow  visualization 


Dye  was  injected  at  the  leading  edge  of  the  delta  wing  just  as  the  freestream 
velocity  is  increased.  Figure  5-2  shows  these  results. 


Figure  5-2  Vortex  breakdown,  dye  visualization 
The  viewing  area  is  from  approxunately  10%  to  60%  of  the  centerline  chord.  The 
vortex  tube  is  viewed  from  the  side,  as  is  indicated  in  the  illustration.  The  dyed  particles 
represent  the  vortex  filaments,  since  the  dyed  particles  are  initially  sucked  into  the  low 
pressure  of  the  cores  of  the  winding  vortex  filaments  and  remain  always  aligned  with  the 

vortex  filaments  by  the  Helmholtz  vortex  theorem  ( ®  =  c  •  8/ ).  The  thick  line  is  the  edge 
of  the  delta  wing.  The  three  images  are  taken  at  successive  time  intervals.  They  are 
chosen  to  show  the  process  of  vortex  breakdown  in  the  vorticity  gradient  region.  In  the 
first  image,  radial  expansion  begins.  In  the  next  image  we  see  continuing  radial 


expansion  and  pile-up.  The  pUe-up  manifests  itself  by  an  increased  concentration  of 
dye.  This  is  why  the  vortex  tube  looks  darker  in  the  radial  expansion  region  in  the  second 
image.  In  the  third  image,  a  turning  point  is  visible.  Backflow  is  also  present,  indicating 
the  presence  of  the  azimuthal  vorticity  sign  switch.  Figure  5-3  shows  the  last  image  m 
more  detail.  All  three  aspects  of  vortex  breakdown  are  visible  in  this  image.  We  see 
radial  expansion  of  the  vortex  tube.  Pile-up  is  present,  as  indicated  by  the  accumulation 
of  dye  near  the  radial  expansion  region.  An  azimuthal  vorticity  sign  switch,  which  will 
be  substantiated  later  by  PIV  data,  is  inducing  backflow.  This  backflow  leads  to  the 
turning  point  and  soon  after,  steady  state  vortex  breakdown. 


(Oe  sign  switch — ^backflow 


Figure  5-3  Dye  flow  visualization,  turning  point 
We  can  compare  this  image  with  a  rendered  image  from  the  vortex  filament 
simulation.  Figure  5-4  shows  this  result.  All  features  of  the  inviscid  self-induction 
mechanism  are  present  in  these  images.  We  see  the  radial  expansion  and  pile-up  along 
the  vortex  tube.  If  we  foUow  the  vortex  filament,  we  see  the  azimuthal  vorticity  sign 
switch.  We  also  see  where  the  vortex  tube  turns  inward  on  itself.  This  is  the  turning 
point.  An  axial  vorticity  sign  switch  characterizes  the  turning  point.  The  axial  vorticity 
sign  switch  causes  coimter-rotation  interior  to  the  vortex  tube  that  makes  the  effects  of 
viscosity  important,  leading  to  turbulences  dowmstream  of  the  turning  point,  as  the 
ejqieriment  shows. 
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Figure  5-4  Vortex  filament  simulation,  features  of  self-induction  mechanism,  T  -  1.5 
Other  dye  flow  visualization  experiments  also  show  the  turning  point.  We  will 
present  these  only  briefly  for  purposes  of  comparison.  These  images  are  from  an 
experiment  with  vortex  breakdown  in  a  pipe  {ref.  8).  This  setup  allows  for  a  much  more 
symmetric  vortex  breakdown. 
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Figure  5-5  Turning  in  vortex  breakdown  experiments  {ref.  8) 


5.1.3  Laser  induced  fluorescence  (LIF)  flow  visualization 

A  fluorescing  dye  is  injected  into  the  vortex  tube  at  the  leading  edge  of  the  delta 
wing.  A  laser  sheet  is  aligned  with  the  vortex  tube  along  its  meridional  plane.  When  the 
laser  sheet  illuminates  the  dye,  the  dye  fluoresces.  Figure  5-6  shows  the  results  of  the 
vorticity  gradient  experiment  using  LIF.  This  is  a  similar  perspective  to  that  of  the  dye 
flow  visualization  except  the  viewing  area  is  smaller.  These  LIF  images  show  the  vortex 
tube  from  40%  to  60%  of  the  centerlme  wing  chord.  Once  again,  the  manifestations  of 
the  vorticity  gradient  are  present.  There  is  radial  ejq)ansion  and  pile-up.  Pile-up  is 
indicated  by  the  increased  intensity  of  the  fluorescence  near  the  radial  expansion  region. 
A  turning  point  is  visible  in  the  second  image.  The  third  image  shows  the  complete 
vortex  breakdown.  The  results  also  support  the  conclusions  of  the  self-induction  theory 
of  vortex  breakdown. 
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Figure  5-6  Laser  induced  fluorescence,  delta  wing  vortex  breakdown 
5.1.4  Particle  image  velocimetry  (PIV) 


Dye  flow  visualization  and  LIF  are  useful  tools  in  visualizing  vortex  breakdown, 
but  they  cannot  tell  us  anything  about  flow  velocity.  Dye  flow  visualization  and  LIF 
show  streaklines,  not  streamlines,  unless  the  flow  is  steady.  However,  we  can  use 
another  tool,  particle  image  velocimetry  (PIV),  to  find  the  velocity  in  a  particular  plane  in 
the  flow.  PIV  works  by  seeding  the  flow  with  particles  and  taking  the  picture  of  these 
particles  at  successive  time  intervals.  Using  a  statistical  correlation  scheme,  PIV 
software  compares  the  images  fi-om  one  to  the  next  and  determines  the  direction  and  how 
far  the  particles  moved.  Since  the  time  between  images  is  known,  we  can  find  the 
velocity  vector  averaged  across  the  time  between  images.  Since  the  time  between  images 
is  usually  very  short,  we  can  say  that  the  velocity  vector  is  almost  an  instantaneous 
picture  of  the  velocity  vectors  in  the  flow  at  that  time.  Like  the  LIF,  a  laser  is  used 
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illuminate  a  plane  within  the  flow.  In  this  case,  the  meridional  plane  of  the  vortex 
tube  is  illuminated. 

The  flow  is  seeded  with  small  hydrogen  bubbles  that  are  lit  up  by  the  laser.  The 
hydrogen  bubbles  are  generated  by  electrolysis  with  electrodes  mounted  at  the  apex  of 
the  delta  wing.  Due  to  their  low  density,  the  hydrogen  bubbles  are  quickly  sucked  into 
the  center  of  the  vortex  tube.  To  mininiize  the  effects  of  buoyancy  of  the  hydrogen 
bubbles,  the  delta  wing  is  mounted  upside  down  within  the  water  timnel.  Also,  because 
the  images  were  taken  at  very  short  time  intervals,  the  effect  of  buoyancy  between  these 
time  intervals  is  negligible.  Figure  5-7  shows  the  setup  of  this  experiment. 


With  this  setup,  we  can  find  the  velocity  vectors  along  in  the  meridional  plane  of 
the  vortex  tube.  We  will  have  axial  and  radial  velocity  data.  With  that  data  we  can  also 
calculate  azimuthal  vorticity. 

We  present  the  PIV  results  at  three  different  times.  The  first  image  shows  the 
vorticity  gradient  region  as  it  first  enters  the  PIV  capture  area.  This  will  be  to.  After  that, 
images  are  taken  every  5  milliseconds.  Three  images  are  then  averaged  together  to  make 
the  effective  time  between  images  1 5  milliseconds.  The  following  results  will  show  the 
images  at  to,  to  +  45  ms,  and  to  +  60  ms. 
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Figvire  5-8  shows  the  velocity  vectors  and  vorticity  contour  plot  at  to.  The 
vectors  in  the  upper  left  show  the  velocity  in  the  meridional  plane  of  the  vortex  tube.  The 
bottom  right  color  contour  plot  shows  the  azimuthal  vorticity.  Here,  the  initiatmg 
vorticity  gradient  is  clearly  visible.  The  darker  red  and  darker  blue  upstream  show  a  clear 


negative  azimuthal  vorticity  gradient 
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We  also  see  radial  expansion  in  the 


velocity  vector  plot,  as  the  vortex  tube  grows  radially.  Pile-up  is  also  present,  since  the 


velocity  decelerates  in  the  axial  direction  along  the  vortex  tube  centerline 
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Figure  5-8  PIV  results,  vorticity  gradient  enters  capture  area,  t  =  to 
Figure  5-9  shows  the  velocity  vectors  and  vorticity  color  contours  at  /  =  /o  +  45 
ms.  Here  the  all  features  of  the  self-induction  mechanism  are  visible.  Radial  expansion 
is  clearly  indicated  in  the  velocity  and  vorticity  color  contour  plots.  Pile-up  is  present 
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and  increasing  in  strength,  as  the  axial  velocity  gradient  becomes  more  negative  in  the 
vorticity  gradient  region.  Now,  we  can  see  the  azimuthal  vorticity  sign  switch,  as 
indicated  in  the  figure.  The  effects  of  this  azimuthal  vorticity  sign  switch  are  evident  in 
the  velocity  vector  plot,  as  the  axial  velocity  decreases  almost  to  zero  in  the  region 
associated  with  the  azimuthal  vorticity  sign  switch. 


Figure  5-9  PIV  results,  azimuthal  vorticity  sign  switch,  r  ==  to  +  45 
Finally,  in  Figure  5-10,  at  r  =  fo  +  60  TWi'  we  see  the  vortex  breakdown  established 
and  beginning  to  move  upstream.  The  azimuthal  vorticity  sign  switch  is  stronger.  The 
maximum  induced  backflow  fi-om  this  azimuthal  vorticity  sign  switch  is  seeking  to 
balance  itself  with  the  relative  freestream  velocity.  The  vortex  breakdown  will  stop 
moving  at  the  location  where  the  induced  backflow  equals  the  relative  fi’eestream 
velocity.  In  the  frame  of  reference  of  these  images  (body-fixed)  this  point  will  be  when 
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the  total  velocity  of  the  particles  is  zero  at  and  around  the  vortex  breakdown  region 
center. 


Figure  5-10  PIV  results,  vortex  breakdown  established,  /  =  to  +  60 
These  PIV  results  capture  the  transient  formation  of  vortex  breakdown,  as  brought 
about  by  a  negative  vorticity  gradient  introduced  in  to  the  vortex  tube  by  an  increase  in 
freestream  velocity.  Here  we  see  all  parts  of  the  self-induction  vortex  breakdown 
mechanism-radial  expansion,  pile-up,  and  azimuthal  vorticity  sign  switching. 

We  can  compare  these  results  directly  with  the  vortex  filament  method  by  looking 
at  the  azimuthal  vorticity  output.  Figure  5-1 1  shows  the  azimuthal  vorticity  results  firom 
T  =  1 .5  in  the  vortex  filament  simulation  at  t  =  +  45  in  the  PIV  experiment.  Many 

features  are  similar,  though  the  experiment  suffers  fi-om  the  lack  of  symmetry  associated 
with  delta  wing  vortex  breakdown.  Still,  the  vorticity  gradient  is  clear,  as  is  the  radial 
expansion  and  azimuthal  vorticity  sign  switch. 
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Azimuthal  Vortictty,  meridonal  plane,  Time:  1 .5 


Figure  5-11  Azimuthal  vorticity,  simulation  (top)  and  experiment  (bottom) 

We  can  also  compare  axial  velocity  profiles  at  different  axial  location.  It  is  well 
known  that  the  axial  velocity  profile  in  a  vortex  tube  is  jet-like  upstream  of  the  vortex 
breakdown.  However,  it  changes  downstream  of  the  vortex  breakdown  to  a  wake-like 
axial  velocity  profile.  If  this  same  change  occurs  in  the  vortex  filament  simulation,  this 
change  in  axial  velocity  fi'om  a  jet-like  to  wake-like  may  be  an  inviscid  mechanism  in 
response  to  the  velocities  induced  by  the  azimuthal  vorticity  in  the  vortex  tube.  Figure 
5-12  compares  the  axial  velocity  at  locations  upstream  and  downstream  of  vortex 
breakdown  with  results  firom  the  PIV  experiment  and  the  vortex  filament  simulation. 
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Figure  5-12  Axial  velocity  profile,  experiment  (top,  middle)  and  simulation  (bottom) 
The  results  are  clearly  similar.  Both  the  vortex  filament  simulation  and  the  PIV 
experiment  show  jet- like  axial  velocity  profiles  upstream  of  the  vortex  breakdown,  and 
wake- like  axial  velocity  profiles  downstream  of  the  vortex  breakdown.  However,  the 
simulation  only  models  inviscid  self-induction.  Thus,  we  can  conclude  that  these 


158 


velocity  profiles  are  caused  by  induced  velocities  and  not  by  viscous  effects.  Figure 
5-13  illustrates  how  this  might  work. 
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Figure  5-13  Azimuthal  vorticity  sign  switch  and  velocity  profile 
If  we  linearly  add  the  freestream  velocity  to  the  velocity  induced  by  azimuthal 
vorticity  components  upstream  of  the  vortex  breakdown,  it  will  enhance  the  freestream 
velocity,  leading  to  the  jet-like  velocity  profile  that  we  see  on  the  left  in  Figure  5-13. 
However,  with  the  azimuthal  vorticity  sign  switch  downstream  of  the  vortex  breakdown, 
the  velocities  induced  by  the  azimuthal  vorticity  will  oppose  the  freestream  velocity, 
leading  to  the  wake-like  velocity  profile  we  see  on  the  right. 

The  dye  flow  visualization,  LIF,  and  PIV  experiments  qualitatively  confirm  the 
features  of  the  self-induction  theory  of  vortex  breakdown.  In  the  experiment,  a  vorticity 
gradient  causes  radial  expansion,  pile-up  and  vorticity  sign  switching,  just  as  in  the 
vortex  filament  simulation.  This  provides  further  evidence  that  vortex  breakdown  is 
initiated  by  a  vorticity  gradient  that  manifests  itself  in  an  inviscid  self-induction  feedback 
mechanism 
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5.2  Case  2,  vortex  breakdown  due  to  vorticUy  gradient  in  vortex  tube  with  fixed 
circulation 

The  author  performed  this  experiment  for  the  purposes  of  investigating  how  vortex 
blowing  could  affect  vortex  breakdown.  The  initial  goal  was  to  delay  vortex  breakdown 
through  periodic  blowing  along  the  leading  edge  of  the  delta  wing  in  order  to  perturb  the 
vortex  core.  However,  it  was  observed  that  this  blowing,  in  fact,  caused  vortex 
breakdown.  The  reason  became  evident  in  light  of  the  self-induction  theory.  If  a  vortex 
tube  is  acted  upon  by  a  pulsed  blowing  that  causes  a  perturbation  along  its  core,  a 
vorticity  gradient  forms  along  the  axis  of  the  core  at  the  perturbation.  Maxwothy,  Mory, 
and  Hopfinger  in  ref.  26  (among  others)  recognize  that  these  types  of  perturbations  in  the 
form  of  waves  can  lead  to  vortex  breakdown,  though  they  do  not  comment  on  the 
importance  of  the  vorticity  gradient  and  self-induction. 

Figure  5-14  illustrates  how  a  perturbation  in  the  vortex  tube  causes  a  negative 
axial  vorticity  gradient  along  the  main  axis  of  the  vortex  tube. 

Perturbed  Vortex  Tube 


H - — r - T - '  ‘ 

0  1  2  3  4 


Axial  location,  z 

Figure  5-14  Vorticity  gradient  caused  by  vortex  tube  perturbation 
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In  this  simplified  example,  the  vorticity  along  the  vortex  tube  is  a  constant 
value,  ®zo,  and  entirely  in  the  axial  direction.  When  the  vortex  tube  is  perturbed,  the 
angle  between  the  tube  and  the  main  axis  of  the  vortex  changes  by  angle  0,  in  the  region 
of  the  perturbation.  This  causes  a  local  change  in  the  axial  vorticity  in  the  region  of  the 
perturbation,  according  to  ’  cos(0) 

This  vorticity  gradient,  if  sufficiently  strong  will  lead  to  vortex  breakdown  in  the 
region  of  that  vorticity  gradient.  In  this  experimental  case,  a  periodic  jet  on  the  leading 
edge  of  the  delta  wing  momentarily  perturbs  the  vortex  on  the  upper  part  of  the  delta 
wing.  Figure  5-15  shows  a  diagram  of  this  setup.  This  perturbation  leads  to  an 
immediate  vortex  breakdown  in  the  region  of  the  perturbation.  A  secondary  effect  (but 
the  main  effect  we  are  showing  here)  is  that  the  lower  vortex  is  also  perturbed.  This  is  a 
weaker  perturbation,  so  the  resulting  vortex  breakdown  happens  slowly  enough  that  the 
breakdown  process  can  be  captured  in  the  following  images.  As  the  perturbation 
convects  downstream,  the  features  described  by  the  self-induction  mechanism  become 
visible.  Figure  5-16  shows  the  dye  flow  visualization  over  time  of  this  perturbation 
leading  to  a  localized  vortex  breakdown. 


Periodic  leading  edge  jet 


Figure  5-16  Vortex  breakdown  in  perturbed  vortex  tube 
This  provides  further  evidence  that  vortex  breakdown  is  initiated  by  a  negative 
vorticity  gradient  and  manifests  itself  in  a  feedback  mechanism  that  includes  radial 
expansion  and  contraction,  pile-up,  and  azimuthal  vorticity  sign  switching  (which  mustt 
be  present  if  there  is  backflow  relative  to  the  freestream  velocity).  These  mechanisms 
culminate  in  the  turning  point  where  there  is  a  sign  switch  in  the  axial  vorticity.  This 
axial  vorticity  sign  switch  is  associated  with  counter  rotating  flows  that  is  very  quickly 
disrupted  due  to  the  formation  of  shear  layers  and  the  resulting  small  scale  turbulence, 
which  is  visible  downstream  of  the  main  breakdown  location. 

5.3  Case  3,  vortex  breakdown  in  cylinder  with  rotating  end  wall 

The  experiment  with  the  delta  wing  in  the  water  tuimel  was  very  useful  for 
showing  the  results  of  an  explicitly  specified  vorticity  gradient.  However,  the  delta  wing 
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disturbs  the  symmetry  of  the  vortex  breakdown,  while  the  vortex  filament  simulation 
is  inherently  symmetric.  We  will  now  look  at  the  results  from  vortex  breakdown  in  a 
cylinder  in  order  to  compare  the  results  with  more  symmetric  cases.  This  ejqjeriment  was 
performed  by  M.P.  Escudier  and  is  published  in  ref.  27.  This  experiment  generates  a 
vortex  by  rotating  the  end  wall  in  a  cylindrical  container  filled  with  water.  This  rotating 
end  wall  causes  rotating  and  recirculation  flow  due  to  Ekman  layer  pumping.  At  the 
center  of  this  rotating  flow,  a  vortex  forms  which  can  undergo  vortex  breakdown.  Figure 
5-17  shows  a  diagram  of  the  experiment  setup. 


Figure  5-17  Vortex  breakdown  in  cylinder,  experiment  setup 
Using  this  setup,  experiments  were  performed  at  various  ratios  of  cylinder  height, 
H,  to  cylinder  radius,  R,  and  at  various  Reynolds  number.  In  this  case,  Reynolds  number 
is  based  on  the  angular  velocity  of  the  rotating  end  wall,  Q,  the  radius  of  the  cylinder,  R, 

Q.jR2 

and  the  kinematic  viscosity  of  the  fluid,  v,  giving  Re  =  — .  Flow  observations  are 

made  using  laser-induced  fluorescence.  The  flow  presented  is  steady,  so  the  streaklines 
from  the  dye  correspond  to  the  streamlines.  Figure  5-18  shows  the  results  at  H/R  =  2  and 


Re  =1854. 


Figxire  5-18  Vortex  breakdown  in  cylinder,  H/R  =  2,  Re  =  1 854 
This  compares  well  with  the  instantaneous  streamline  plot  from  the  present  vortex 
filament  method,  as  presented  in  Figure  5-20.  In  both  results,  the  radial  expansion,  pile- 
up,  and  backflow  near  the  main  vortex  breakdown  location  is  evident.  Also,  what  is 
apparent  is  a  second  region  of  radial  expansion  downstream  of  the  main  breakdown 
location.  This  is  clear  in  Escudier’s  experiment  and  also  in  the  vortex  filament 
simulation  result  as  Figure  5-20  illustrates.  This  indicates  the  formation  of  a  second 
vortex  breakdown  downstream  of  the  first.  Since  the  vortex  filament  simulation  captures 
this  result,  it  too,  is  likely  caused  by  a  vorticity  gradient,  working  through  the  inviscid 
feedback  mechanism.  This  can  be  described  by  the  self-induction  theory  of  vortex 
breakdown.  Due  to  another  ©e  sign  switch  associated  with  the  radial  contraction  of  the 
vortex  tube,  as  seen  in  Figure  5-19,  the  sign  of  ©e  returns  to  its  initial  value,  and  through 
the  repetitive  inviscid  mechanism,  the  secondary  vortex  breakdown  and  the  secondary  ©e 
appears.  The  capturing  of  the  secondary  vortex  breakdown  by  the  inviscid  simulation 
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may  imply  that  the  inviscid  process  can  persist  downstream  of  the  primary  vortex 
breakdown. 


Figure  5-19  Azimuthal  vorticity  and  secondary  vortex  breakdown,  T  -  2.4 


Streamlines,  meridional  plane,  Time;  2.4 


Figure  5-20  Instantaneous  streamlines,  vortex  filament  simulation,  T  -  2.4 

5.4  Case  4,  NavierStokes  solution  to  case  2 

The  Escudier  e)q)eriment  in  the  cylinder  provides  only  streakline  results. 
Masanori  Kikuchi  (ref.  28)  took  the  conditions  of  this  experiment  and  calculated  a 
numerical  solution  using  the  full  Navier-Stokes  equation.  The  benefit  of  this  approach  is 


165 


that  we  can  get  data  that  the  physical  experiment  cannot  provide.  Because  of  the 
relatively  low  Reynolds  numbers  involved,  the  Navier-Stokes  equations  can  be  solved  to 
a  very  precise  degree  of  accuracy,  giving  us  confidence  in  the  solution.  The  Navier- 
Stokes  equations  capture  all  possible  physics  in  this  problem,  including  the  effects  of 
viscosity.  However,  if  the  results  from  the  Navier-Stokes  solution  compare  well  with  the 
vortex  filament  simulation,  which  does  not  include  the  effects  of  viscosity,  it  can  be 
concluded  that  the  vortex  breakdown  process  is  nearly  inviscid. 

The  results  for  the  numerical  Navier-Stokes  solution  for  the  Escudier  e3q)eriment 
are  presented  in  the  following  figures.  These  are  for  the  case  of  H/R  =  2.5  and  Re  = 
2765.  Figure  5-21  shows  the  streamlines  from  the  numerical  solution  to  the  Navier- 
Stokes  Equation  (L)  and  Escudier’s  experiment  (R).  The  results  compare  well  with  each 
other,  as  one  would  expect.  Radial  expansion  and  contraction  and  backflow  are  clear  in 
both  results.  The  turning  point  is  clearly  visible  in  Escudier’s  experimental  results, 
though  it  is  not  as  clear  in  the  simulation  results.  Still,  secondary  breakdown  exists  in 
both  cases.  Overall,  the  Navier-Stokes  simulation  compares  well  with  the  results  of  the 
experiment.  Thus,  we  can  be  confident  that  the  data  calculated  with  the  numerical 
solution  will  also  be  accurate. 
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Figure  5-21  Streamlines,  Navier-Stokes  Simulation  (L)  and  Experiment  (R) 

Figure  5-22  shows  the  azimuthal,  axial,  and  radial  vorticity  data  for  the  Navier- 
Stokes  simulation  of  the  Escudier  experiment.  These  results  present  the  data  on  the 
meridional  half-plane.  Here  we  see  many  of  the  same  featmes  as  we  did  in  the  vortex 
filament  simulations.  Namely,  vorticity  sign  switching  is  clearly  present  in  the  azimuthal 
and  radial  vorticity.  In  the  axial  vorticity,  there  is  no  sign  switching,  though  there  is  a 
marked  decrease  in  axial  vorticity  in  the  vortex  breakdown  regioa  This  may  be  due  to 
the  weak  axial  flow  rate  and  the  constraints  of  the  cylinder  walls  in  close  proximity  to  the 
vortex  breakdown  location.  It  may  also  be  that,  with  viscous  dissipation  possible,  any 
axial  vorticity  sign  switch  is  quickly  worn  away  by  these  viscous  forces,  as  was 
mentioned  in  the  original  discussion  of  the  turning  point.  Even  if  this  is  the  case,  this 
viscous  action  does  not  cause  vortex  breakdown;  it  simply  wears  away  some  of  the 
indicators  of  vortex  breakdown,  after  that  breakdown  has  occurred. 
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Figure  5-22  Meridional  Plane  Data,  Navier-Stokes  Simulation,  a)e  (L),  ©z  (C),  cOr  (R) 

From  the  Navier-Stokes  solution,  a  three-dimensional  plot  of  vortex  lines  can  be 
constructed,  as  shown  in  Figure  5-23.  This  is  different  than  the  way  vortex  lines  are 
constructed  in  the  vortex  filament  simulation,  since  in  the  vortex  lines  are  implicit  the 
simulation  in  the  vortex  filament  method,  while  in  a  Navier-Stokes  solution  they  are 
foimd  through  post-processing  of  the  raw  data.  Still,  the  results  compare  well  with  the 
vortex  filament  simulation  in  terms  of  the  self-induction  theory  of  vortex  breakdown.  In 
the  Figure  5-23  result,  pile-up  and  radial  expansion  are  clear.  Also,  the  location  of  the 
azimuthal  vorticity  sign  switch  is  visible  when  you  decompose  the  vortex  line  into  vector 
components.  Compare  the  result  firom  the  Navier  Stokes  numerical  solution  to  the  vortex 
filament  method  results  at  T  =  0.9,  as  in  Figure  5-24. 


©e  sign  switch 
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While  there  is  no  axial  vorticity  sign  switch  in  the  Navier-Stokes  solution  to 
Escudier’s  experiment,  there  are  still  similarities  between  the  axial  vorticity  profile  in  the 
Navier-Stokes  simulation  and  the  vortex  filament  method  simulation.  Figure  5-25 
presents  the  axial  vorticity  along  the  vortex  tube’s  axis  at  various  radial  locations.  The 
trends  in  axial  vorticity  are  very  similar,  with  the  peak  axial  vorticity  occurring  along  the 
centerline.  That  peaks  in  axial  vorticity  decreases  near  the  vortex  breakdown  region 
before  increasing  and  then  retirniing  to  zero  downstream  of  the  breakdown  location.  The 
peaks  decrease  farther  away  fi-om  the  centerline  in  both  simulations.  The  jaggedness  in 
the  vortex  filament  simulation  solution.  Figure  5-25  (R),  is  due  to  the  resolution  of  the 
grid  that  the  vorticity  data  is  calculated  on.  Calculating  this  data  is  not  implicit  to  the 
vortex  filament  method,  though  it  is  for  the  Navier-Stokes  simulation.  Still,  the  trends 
compare  well  with  each  other,  with  the  major  difference  being  the  absence  of  the  axial 
vorticity  sign  switch,  which  may  be  the  result  of  viscosity. 
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Figure  5-25  Axial  vorticity  at  various  radial  location  along  vortex  tube,  Navier-Stokes 
Simulation  (L)  and  vortex  filament  simulation  (R) 

5.5  Other  e^erimental  and  computation  evidence 

Other  experiments  and  computational  studies  give  support  to  the  components  of 
the  self-induction  theory  of  vortex  breakdown.  They  provide  specific  examples  of  pile- 
up,  radial  expansion,  and  azimuthal  and  axial  vorticity  sign  switching. 

Experiments  with  vortex  flows  and  shock  waves  show  that  vortex  breakdown 
occurs  immediately  downstream  of  the  location  where  a  vortex  crosses  a  shock  wave  (see 
ref.  29  for  details).  Figure  5-26  shows  a  diagram  of  this.  In  terms  of  the  self-induction 
theory  of  vortex  breakdown,  the  shock  causes  a  strong  axial  vorticity  gradient  in  the 
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vortex  tube,  as  a  result  of  the  axial  velocity  decrease  downstream  of  the  shock.  This 
also  causes  pile-up  within  the  vortex  tube.  All  of  this  begins  the  inviscid  feedback 
mechanism  that  leads  to  steady  state  vortex  breakdown. 

Rtelreutatlon  Subtilt  B^cotnpr^Mon 


Figure  5-26  Shock  wave/vortex  interactions  (ref.  29) 

A  computational  study  by  M.R.Visbal  provides  detailed  vorticity  data  for  vortex 
breakdown  (see  ref  30  for  details).  He  compares  his  computational  results  to 
experiments  by  Shih  and  Ding.  These  results  shows  the  azimuthal  vorticity  sign  switch 
in  the  region  of  vortex  breakdown,  confirming  the  same  results  in  the  vorticity  gradient 
experiment  of  S.  Srigrarom  (ref  4)  and  the  current  vortex  filament  simulation.  Figure 
5-27  shows  these  results. 
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Figure  5-27  Azimuthal  vorticity  sign  switch  {ref.  30) 

Visbal’s  computational  results  also  show  the  axial  vorticity  sign  switch.  Figure 
5-28  shows  these  results  in  both  the  meridional  and  cross  flow  plane.  Visbal  interprets 
these  results  by  saying,  “ . . .  negative  axial  vorticity  concentrations  may  be  the  result  of 
the  severe  upstream  tilting  of  part  of  the  incoming  vorticity  by  the  strong  reverse  flow 
induced  by  the  spiral”  (ref  30).  This  is  exactly  what  the  vortex  filament  simulation 


shows. 


Figure  5-28  Axial  vorticity  sign  switch,  meridional  (top)  and  cross  flow  (bottom)  plane 

{ref.  30) 

The  experimental  data  and  the  results  from  other  computer  simulations  compare 
well  with  the  results  from  the  vortex  filament  method.  These  results  also  fit  well  into  the 
self-induction  theory  of  vortex  breakdown.  This  gives  further  support  to  that  theory, 
indicating  that  vortex  breakdown  is  an  inviscid  phenomenon,  caused  by  a  vorticity 
gradient,  and  brought  about  by  an  inviscid  feedback  mechanism  consisting  of  radial 
expansion  and  contraction,  pile-up,  and  vorticity  sign  switching.  With  this  knowledge  of 
the  mechanism  that  may  lead  to  vortex  breakdown,  one  can  better  formulate  ways  to 
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weaken,  prevent,  or  delay  vortex  breakdown  and  its  corresponding  negative  effects  on 
aircraft  flight. 

5.6  Vortex  breakdown  control  experiments  and  the  self-induction  theory 

The  self-induction  theory  of  vortex  breakdown  provides  a  clearer  picture  of  the 
unsteady  formation  of  vortex  breakdown.  It  describes  an  inviscid  feedback  mechanism 
that  starts  fi-om  a  negative  vorticity  gradient  and  includes  radial  expansion  and 
contraction,  pile-up,  and  vorticity  sign  switching.  If  a  negative  vorticity  gradient  and 
these  mechanisms  lead  to  vortex  breakdown,  then  the  interruption  of  one  or  more  of  these 
mechanisms  may  delay,  weaken,  or  prevent  vortex  breakdown.  Thus  we  can  control 
vortex  breakdown  by  affecting  the  vorticity  gradient,  the  radial  expansion  and 
contraction,  pile-up,  or  the  vorticity  sign  switching.  Various  researchers  have 
successfully  manipulated  the  location  of  vortex  breakdown  over  delta  wings  and  other 
experimental  apparatus.  We  will  look  at  several  of  these  e>q)eriments  and  discuss  them  in 
terms  of  the  self-induction  theory  of  vortex  breakdown. 

5. 6. 1  Vorticity  gradient 

If  vortex  breakdown  is  initiated  by  a  vorticity  gradient,  then  anything  that  would 
change  that  vorticity  gradient  should  change  the  location  of  the  vortex  breakdown.  As 
was  previously  described,  a  velocity  gradient  over  a  wing  can  lead  to  a  vorticity  gradient. 
Thus,  if  that  velocity  gradient  (or  the  associated  adverse  pressure  gradient)  is  modified, 
the  vorticity  gradient  will  also  be  modified.  A  weaker  vorticity  gradient  leads  to  a  more 
downstream  vortex  breakdown  location,  so  a  weaker  velocity  gradient  should  delay 
vortex  breakdown. 

This  is  what  happens  in  the  experiments  of  A.  Mitchell  et  al.  See  ref.  31  for 
details.  Through  symmetrical  and  asymmetrical  blowing  along  the  trailing  edge  of  a 
slender  delta  wing,  they  were  able  to  delay  vortex  breakdown.  In  terms  of  the  self- 
induction  theory,  this  can  be  understood  in  terms  of  a  change  in  the  vorticity  gradient. 
Trailing  edge  blowing  changes  the  velocity  profile  over  the  delta  wing  and  weakens  the 
adverse  pressure  gradient  over  the  leeward  surface  of  the  wing.  This  corresponds  to  a 
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decrease  in  the  velocity  gradient  and  a  corresponding  decrease  in  the  vorticity 
gradient.  As  the  self-induction  theory  predicts,  the  location  of  vortex  breakdown,  at  a 
given  angle  of  attack,  should  move  downstream  when  the  trailing  edge  blowing  is  active. 

Another  way  to  change  the  vorticity  gradient  is  to  change  the  shape  of  the  leading 
edge  near  the  leading  edge  tip.  This  is  where  most  circulation  in  the  vortex  tube 
originates,  so  it  will  be  the  most  sensitive  to  small  changes  in  configuration.  Adding 
fillets  or  rounding  off  sharp  edges  will  reduce  the  circulation  in  the  vortex  tube,  while 
leaving  larger  flow  properties  like  the  velocity  gradient  over  the  wing  and  the  strength  of 
the  vortex  sheet  unchanged.  Thus,  the  total  amount  of  circulation  is  reduced,  and  the 
maximum  possible  induced  backflow  velocity  is  also  decreased,  resulting  in  a  vortex 
breakdown  location  either  further  downstream  over  a  wing  at  a  fixed  angle  of  attack,  or 
increasing  the  critical  angle  of  attack  where  vortex  breakdown  will  first  occur  at  a  given 
location  over  the  wing. 

This  is  what  is  done  by  S.  Srigrarom  in  an  experiment  designed  to  use  the  results 
of  the  self-induction  theory  to  delay  vortex  breakdown.  See  refA  for  details  about  this 
experiment.  The  experiment  compares  the  breakdown  location  over  an  F/ A- 18  model  (a) 
with  the  leading  edge  of  the  strake  left  sharp  and  (b)  with  the  leading  edge  rounded  using 
clay.  Figure  5-29  shows  the  modification  to  the  F/A-18  model.  This  rounded  leading 
edge  should  reduce  the  circulation  in  the  vortex  tube  and  therefore  delay  vortex 


breakdown. 
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Left  side,  NO  clay  attached  Right  side,  clay  attached  at  the 

juncture  of  LEX  and  the  body. 

Figure  5-29  F/A-18  Leading  edge  strake  modification  (ref.  4) 

The  modified  F/A-18  is  placed  in  a  water  tunnel  at  25  degrees  angle  of  attack  with 
a  tunnel  speed  of  5  cm/s.  The  aircraft  is  aligned  such  that  there  is  no  yaw  angle.  Dye  is 
used  for  flow  visualization.  The  vortex  breakdown  location  is  noted  by  rapid  radial 
expansion,  backflow,  and  downstream  turbulence.  Figure  5-30  shows  the  results.  The 
image  on  the  left,  with  the  sharp  leading  edge  extension,  shows  the  normal  vortex 
breakdown  location  for  this  angle  of  attack  at  approximately  0.3c.  However,  with  the 
rounded  leading  edge  (right  image)  there  is  no  visible  vortex  breakdown.  There  is  only  a 
disruption  of  the  vortex  tube  when  it  encounters  the  vertical  tail. 


177 


Clean  strake  LE  at  root  Rounded  strake  LE  at  root 


Figure  5-30  Delay  in  vortex  breakdown  location  with  roimded  strake  LE  {ref.  4) 
Reducing  the  circulation  in  the  vortex  tube  may  delay  vortex  breakdown,  but  it 
will  also  reduce  the  high  axial  velocities  and  corresponding  low  pressures  that  occur 
upstream  of  the  vortex  breakdown  location  and  provide  the  increase  in  lift  associated 
with  vortex  flows  over  wings.  Still,  the  vortex  will  exist  intact  over  a  larger  area,  and  the 
downstream  turbulence  that  leads  to  high-cycle  fatigue  will  be  reduced,  so  the  overall 
effect  may  still  be  positive.  The  net  effect  will  require  further  testing. 

5. 6. 2  Radial  Expansion 

Radial  expansion  is  one  of  the  first  manifestations  of  a  vorticity  gradient  in  a 
vortex  flow.  If  this  radial  expansion  can  be  limited  or  delayed,  it  will  delay  the  formation 
of  vortex  breakdown.  In  a  vortex  generated  in  a  tube,  the  walls  of  the  tube  constrain  the 
radial  expansion.  If  the  radius  of  the  tube  is  sufficiently  small,  as  compared  to  the  vortex 
core  diameter,  radial  expansion,  and  in  turn,  vortex  breakdown  can  be  delayed. 

Similarly,  over  a  delta  wing,  a  vortex  sheet  develops  that  constrains  the  radial  expansion 
of  the  vortex  core.  Figure  5-31  shows  this  vortex  sheet  as  it  separates  fi-omthe  leading 
edge  of  a  delta  wing. 
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Figure  5-31  Vortex  flow  over  a  delta  wing  (ref.  9) 

The  difference  between  this  vortex  sheet  and  the  walls  of  a  pipe  is  that  the  vortex 
sheet  is  not  solid.  Thus  the  strength  of  the  vortex  sheet  will  determine  what  amount  of 
radial  expansion  can  occur  in  the  vortex  core.  If  a  given  vortex  sheet  with  a  certain 
strength  allows  for  sufficient  radial  expansion  for  vortex  breakdown  to  occur  at  one 
location  over  a  delta  wing,  an  increase  in  strength  of  that  vortex  sheet  can  delay  the 
location  of  vortex  breakdown  by  reducing  the  radial  expansion  mechanism.  Reducing 
this  radial  expansion  reduces  pile-up.  Together,  they  reduce  the  vorticity  sign  switching, 
and  the  maximum  possible  induced  backflow  velocity,  thereby  delaying  vortex 
breakdown. 

This  seems  to  be  the  case  in  the  experiment  by  Guy,  Morrow,  and  McLaughlin. 
Through  periodic  blowing  along  the  leading  edge  of  a  delta  wing,  they  delay  vortex 
breakdown  by  up  to  0.2c  (ref  32).  The  effect  of  periodic  blowing  is  mostly  jet-like,  as 
the  suction  stroke  tends  to  draw  from  low  energy  flow  in  the  boimdary  layer  near  the 
blowing  ports.  Therefore,  periodic  blowing  along  the  leading  edge  acts  to  strengthen  the 
vortex  sheet.  This  strengthened  vortex  sheet  constrains  radial  expansion  and  delays 
vortex  breakdown. 
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One  might  wonder,  while  increasing  angles  of  attack  should  also  strengthen 
the  vortex  sheet  and  thereby  delay  vortex  breakdown,  why  breakdown  tends  to  occur 
further  upstream  with  increased  angle  of  attack.  The  answer  lies  in  understanding  that,  as 
angle  of  attack  increase,  circulation  in  the  vortex  sheet,  and  thereby  the  vortex  tube,  also 
increases.  However,  the  adverse  pressure  gradient  and  the  associated  the  negative  axial 
velocity  gradient  along  the  chord  of  the  wing  also  increases,  thereby  increasing  the 
strength  of  the  negative  vorticity  gradient.  Thus,  the  effects  of  the  increasing  the  angle  of 
attack  is  to  not  only  impose  a  stronger  radial  constraint,  but  also  to  strengthen  the 
negative  vorticity  gradient  and  the  induced  velocities  that  cause  radial  expansion. 
However,  this  vorticity  gradient  also  causes  pile-up  and  vorticity  sign  switching.  Thus, 
the  overall  effect  of  increasing  angle  of  attack,  while  tending  to  delay  vortex  breakdown 
in  terms  of  radial  expansion,  is  to  move  the  vortex  breakdown  location  upstream  in  terms 
of  pile-up  and  vorticity  sign  switching.  However,  when  these  mechanisms  can  be 
decoupled,  like  when  the  strength  of  the  radially  constraining  vortex  sheet  increases 
independently  of  vorticity  gradient  by  means  such  as  leading  edge  blowing,  vortex 
breakdown  can  be  delayed.  Aside  from  weakening  the  negative  vorticity  gradient, 
decoupling  the  inviscid  feedback  mechanism  of  radial  expansion,  pile-up,  and  vorticity 
sign  switching,  will  be  an  effective  way  of  weakening,  delaying,  or  preventing  vortex 
breakdown. 

5.6.3  Pile-up 

The  effects  of  pile-up  can  be  delayed  by  adding  an  axial  velocity  component 
along  the  axis  of  the  vortex  tube.  Doing  this  will  implicitly  have  the  effect  of  delaying 
vortex  breakdown  by  requiring  an  increased  maximum  induced  backflow  velocity  to 
equal  the  now  increased  downstream  axial  velocity.  This  may  require  a  more  developed 
(and  hence  further  downstream)  vortex,  or  may  not  be  possible  for  a  given  amount  of 
circulation  in  the  vortex  tube.  Still,  axial  blowing  can  be  interpreted  in  terms  of 
suppressing  pile-up,  thereby  weakening  radial  expansion  and  vorticity  sign  switching. 
The  overall  result  is  to  move  the  location  of  vortex  breakdown  over  a  delta  wing 
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downstream  or  to  delay  the  formation  of  vortex  breakdown  for  a  given  chordwise 
wing  location  until  higher  angles  of  attack.  In  their  experiments  with  axial  blowing  over 
a  delta  wing  A.  Mitchell  et  al  showed  that  they  could  delay  vortex  breakdown  {ref.  33). 
Increases  in  blowing  along  the  core  ftirther  delay  vortex  breakdown,  as  the  self-induction 
theory  would  predict. 

5. 6. 4  Vorticity  sign  switching 

The  azimuthal  vorticity  sign  switching  is  a  major  factor  in  the  development  of  vortex 
breakdown  because  it  allows  for  induced  backflow.  This  backflow  leads  to  feedback 
with  pile-up  and  radial  expansion.  Vorticity  sign  switching  occurs  when  a  portion  of  the 
vortex  tube  rotates  more  slowly  than  the  rest  of  the  vortex  tube.  This  causes  the  vortex 
filaments  to  tilt,  with  the  azimuthal  component  of  the  vorticity  vector  changing  sign. 

This  could  be  delayed  or  weakened  if  an  addition  amount  of  azimuthal  velocity  is  added 
to  the  vortex  tube.  Conversely,  breakdown  could  be  enhanced  by  somehow  slowing 
down  the  rotation  of  the  vortex  tube  in  the  region  of  the  vorticity  gradient.  Experiments 
by  H.  Husain  et  al  show  just  this  result.  See  ref  34  for  details.  In  their  experiments,  a 
vortex  is  generated  in  a  cylinder  in  a  setup  similar  to  the  previously  described  Escudier 
experiment.  However,  a  thin  rod  is  added  along  the  centerline  of  the  cylinder  where  the 
vortex  tube  that  will  develop  due  to  the  rotation  of  the  end  wall.  That  rod  could  be 
rotated  in  the  same  direction  as  the  vortex,  or  it  could  be  counter-rotated  with  respect  to 
the  main  rotating  flow.  With  the  rod  stationary,  vortex  breakdown  is  observed  as  usual. 
When  the  rod  is  co-rotated  at  a  sufficient  rate,  vortex  breakdown  is  no  longer  be  observed 
in  the  flow.  However,  when  the  rod  is  counter-rotated  with  respect  to  the  main  vortex 
flow,  the  vortex  breakdown  was  enhanced,  with  increased  radial  expansion  and 
strengthened  recirculation.  This  can  be  interpreted  in  terms  of  the  azimuthal  vorticity 
sign  switch.  The  co-rotating  rod  added  a  rotation  velocity  to  the  vortex  core  large  enough 
to  overcome  the  rotation  rate  difference  initiated  by  radial  expansion.  Thus,  since  this 
change  in  rotation  rate  was  decreased,  the  azimuthal  vorticity  sign  switch  could  not 
develop  to  a  sufficient  strength  to  cause  vortex  breakdown.  On  the  other  hand,  when  the 


rod  counter-rotated,  it  enhanced  the  azimuthal  sign  switching,  leading  to  increased 
backflow  and  radial  expansion.  Other  ways  to  exploit  this  effect  might  include  blowing 
tangential  to  the  vortex  tube  or  wing  surface  to  increase  the  azimuthal  velocity  and 
weaken  or  remove  the  possibility  of  azimuthal  vorticity  sign  switching.  Zhi  Shi  did  just 
this  to  delay  vortex  breakdown  over  a  delta  wing.  See  ref.  35  for  details.  Figure  5-32 
shows  the  setup  of  the  jets  on  the  surface  of  a  delta  wing. 


Figure  5-32  Tangential  blowing  setup  (ref  35) 
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In  this  setup,  jet  three  blows  tangential  to  the  vortex  tube  in  the  azimuthal 
direction.  With  this  blowing  turned  on,  the  azimuthal  vorticity  sign  switching,  and  hence 
vortex  breakdown  is  delayed.  Figure  5-33  shows  these  results.  The  variable  is  the  jet 
momentum  coefficient,  X/C  is  the  normalized  chordwise  breakdown  location,  and  a  is 
angle  of  attack. 


X/C 


a=17.5-^ 

a-=25^ 

0=30*^ 


Cm 


The  variation  of  breakdown  position  with  tangential 
jet  blowing  (model  II,  jet  3,  Vqo  =  2.73  in/sec). 


Figure  5-33  Delaying  vortex  breakdown  with  tangential  jet  (ref.  35) 


183 


6.  Suggested  future  work 

Much  work  remains  to  be  done  on  this  subject.  Specifically,  we  need  to  work  to 
quantify  these  vortex  filament  method  results  and  look  for  the  limiting  values  of  vorticity 
gradient  that  will  lead  to  vortex  breakdown.  Furthermore,  the  case  with  constant 
circulation  and  a  velocity  gradient  needs  to  be  explored.  This  is  likely  the  most  physical 
way  of  developing  a  vorticity  gradient.  One  way  to  do  this  would  be  to  treat  the 
freestream  velocity  as  a  function  of  axial  distance,  z,  and  time,  t.  In  further  research,  we 
need  to  apply  the  various  possible  techniques  for  delaying  vortex  breakdown  to  the 
vortex  filament  simulation.  To  aid  in  doing  this,  the  simulation  program  could  be  written 
in  a  generalized,  user-fi-iendly  interface  that  would  allow  it  to  become  a  type  of  vortex 
filament  “lab”  where  various  initial  conditions,  boundary  conditions,  and  control  devices 
could  be  setup  easily  to  determining  their  effect  on  the  resulting  vortex  motion. 
Ultimately,  it  is  the  hope  that  this  research  will  lead  to  better  methods  for  controlling 
vortex  breakdown  for  the  purposes  of  flow  and  mixing  control. 
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7.  Conclusions 

Vortex  flows  are  common  in  nature  and  in  aircraft  flows.  The  benefit  of  vortex 
flows  generated  by  strakes  or  swept  wings  on  aircraft  is  their  ability  to  energize  the  flow, 
thereby  delaying  stall  and  increasing  the  maximum  lift  coefficient.  However,  vortex 
flows  can  undergo  rapid  radial  expansion  accompanied  by  backflow  and  turbulence  that 
has  decidedly  negative  effects  on  an  aircraft,  causing  instabilities  when  the  vortex 
breakdown  location  varies  asymmetrically,  and  structural  problems  in  terms  of  high  cycle 
fatigue — a  problem  especially  prevalent  on  twin-tailed  aircraft  like  the  F/A-18.  While 
there  has  been  much  research  aimed  at  understanding  vortex  breakdown,  a  clear,  concise, 
and  widely  accepted  theory  is  still  not  available.  Furthermore,  the  theories  that  are 
available,  while  useful  m  describing  flows  that  eiqperience  vortex  breakdown,  are  not 
very  instructive  in  terms  of  specifying  ways  to  control  vortex  breakdown.  It  is  the 
purpose  of  this  research  to  provide  of  framework  for  imderstanding  the  initial  cause  of 
vortex  breakdown  that  will  lead  researchers  to  effective  methods  to  control  that 
breakdown.  The  self-induction  theory  of  vortex  breakdown  provides  this  framework. 

The  self-induction  theory  of  vortex  breakdown  is  studied  using  computer 
simulations  and  experiments,  as  well  as  examining  the  results  of  past  experiments.  The 
vortex  filament  simulation  and  the  vorticity  gradient  experiments,  in  additions  to  the 
experiments  in  vortex  breakdown  control  provide  strong  evidence  that  vortex  breakdown 
is  an  inviscid  mechanism,  initiated  by  a  negative  vorticity  gradient.  Just  as  in  the  vortex 
filament  simulation,  these  experiments  show  that  a  vorticity  gradient  initiates  the  process 
of  radial  expansion,  pile-up,  and  vorticity  sign  switching  that  leads  to  vortex  breakdown. 
This  provides  fiuther  evidence  that  vortex  breakdown  is  initiated  by  a  vorticity  gradient 
that  manifests  itself  as  a  radial  change  in  flow  through  an  inviscid  self-induction  feedback 
mechanism.  This  self-induction  theory  shows  how  a  feedback  mechanism  works  through 
radial  expansion,  pile  up,  and  vorticity  sign  switching,  to  form  a  turning  point,  thereby 
producing  all  the  features  of  vortex  breakdown. 
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The  primary  benefit  of  the  self-induction  theory  of  vortex  breakdown  is  that  it 
allows  for  a  clearer  understanding  of  mechanism  that  mitiates  vortex  breakdown  while 
also  presenting  clear  ways  of  preventing,  weakening,  or  delaying  that  breakdown. 
Previous  theories  of  vortex  breakdown  do  not  illuminate  possible  control  solution  as 
much  as  they  identify  the  existence  of  vortex  breakdown.  With  this  new  paradigm  of  the 
self-induction  of  vortex  breakdown,  researchers  have  clear  guidelines  to  supplement  their 
engineering  intuition  about  how  to  control  vortex  breakdown.  These  include  changing 
the  vorticity  gradient,  applying  or  removing  radial  constraints,  increasing  or  decreasing 
pile-up,  and  preventing  or  enhancing  the  azimuthal  vorticity  sign  switching.  The 
experiments  that  have  been  presented  provide  evidence  for  this  way  of  thinking  of  vortex 
breakdown.  These  past  successes  in  controlling  vortex  breakdown  point  to  the  self- 
induction  theory  of  vortex  breakdown,  which  in  turn  provide  a  clear  framework  for  future 
successes  in  manipulation  and  control  of  vortex  breakdown. 
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8.  APPENDIX  A — Vortex  Filament  Simulation  Source  Code 


8.1  vlinec.cpp  -  Main  simulation  code 

1 1 - 

#include  <vcl.h> 
tinclude  <f stream. h> 

#include  <istream.h> 

#include  <iostream.h> 

#include  <conio.h> 

#include  <fastmath.h> 

#include  <time.h> 


#pragma  hdrstop 

//initialize  subroutines 
int  read_data  ( )  ; 
int  init  () ; 
void  cond_out (char*) ; 

int  out (double**,  double**,  double**,  char*,  int,  int); 

void  calc_ind  ( ) ; 

void  r_calc_ind  ( ) ; 

void  init_ind ( ) ; 

void  u_update ( ) ; 

void  sym_calc ( ) ; 

void  inir|lncasrrad_spoke( int, double**,  double**,  double**, double**, double, 
int , int* , int , int , double* , double , double , double , double , double , 
double, double, double** , double**, double** , double**, int, int, int) ; 


II - Declare  Global  Variables - - 

int  in  case, velvor_on, blow_on, n_f il, wall_on, wall_ax__node, wall_circ__node, 
double’*t_final, U_inf ,beta, r_tubeO, r_tubel, r_tube2, r_tube3, r_tube4, r_tube5, 
r_tube6,  r_tube7,  r_tube8,  r_tube9,  r_tubelO,  r_tubell,  r_tubel2; 
double  sigma_init,C, alpha, U; 

double  z  on, z_init,x_init,y_init,dz,dt,t_v_out,x_cent, z_start; 

int  t  index,  *on__node,  max_node, out_s kip,  t_out_s tart,  out__count, x_nodes, 

z  nodes ,  end_f ix ,  * r_f il ,  n_rad,  r__f il^in,  r_n_f il ,  r_max_node ,  r_on_f il ; 
double  t_real,pi, *core_con, *r_tube; 

int  n_spokes;  ^  .  +. 

char  *r_node_out,  *r_sigma_out,  *gainma_out,  *r_gamma_out,  *node_out,  *signia__out; 
double  u  induced_x,u_induced_y,u_induced_z; 

double  **x  n,  **y_n,  **z_n,  **u_induced_total_x,  **u_induced_total_y, 

**u_induced_total_z,  **signia; 

double  *  *  gamma ,  start_C_inc ,  s  top_C_inc ,  C_min ; 

double  **r  gamma, **r_x_n,  **r_y_n,  **r_z_n,  **r_u_induced_total_x, 
**r_u_induced_total_y,  **r_u_induced_total_z,  **r_sigma, **r_core_con, 
**r_seg_length; 

int  main (int  argc,  char  **argv) 

{ 

//define  output  files 
node_out  =  "vlinec.dat"; 
sigma_out  =  "sigma__out.dat”; 
r_node_out  =  "r_vlinec.dat”; 
r_sigma_out  ==  "r_sigma_out.dat"; 
gamma__out  =  "gamma_out.dat"; 
r  gamma_out  =  "r_gamma_out.dat”; 


int  cl  =  clock  0; 

read  dataO;  //read  input  file 


initO ///initialize  program  variables 


//output  to  screem 

cout  «  "USE  FOR  CASES  THAT  ARE  ALWAYS  SYMMETRIC  ONLY!"  «endl; 

cout  «  "T_final:  "«t_final  «  endl; 

cout  «  "Total  iterations:  "«t_final/dt  «  endl; 

cout  «  "Max  nodes:  "  «max_node  «  endl; 

cout  «  "U_inf:  "«U_inf  «  endl; 

cout  «  "dz:  "«d2  «  endl; 

cout  «  "Initial  condition  case:  "«in_case  «  endl; 

//initialize  output  files,  will  overwrite  previous  files! 
ofstream  out file (node_out,  ios::trunc); 
if  (loutfile)  { 

cout  «"Error  opening  output  file"«endl;} 
outf ile . close  ( ) ; 

ofstream  sigmafile (sigma_out,  ios::trunc); 
if  (! sigmafile)  { 

cout  «”Error  opening  output  file"«endl;} 
sigmafile . close ( ) ; 

ofstream  r_outfile (r_node_out,  ios::trunc); 
if  (!r_outfiie)  { 

cout  «"Error  opening  output  file"«endl;} 
r_outf ile . close { ) ; 

ofstream  r_sigmafile (r_sigma_out,  ios::trunc); 
if  ( !r_sigmafile)  { 

cout  «"Error  opening  output  file"«endl;} 
r_sigmaf ile . close ( ) ; 

ofstream  gamma file (gamma_out,  ios : : trunc) ; 
if  (igammafile)  { 

cout  «"Error  opening  output  file"«endl;} 
gammaf ile . close  ( ) ; 

ofstream  r_gammafile (r_gamma__out,  ios: : trunc); 
if  ( ! r_gammafile)  { 

cout  «”Error  opening  output  file"«endl;} 
r_gammaf ile . close ( ) ; 

//output  simulation  time  to  screen 

cout  «  "t_index:  "«  t  index  «  "  "  «  '■t_real:  "«t_real  «  endl; 

//write  program  inputs  to  1st  row  of  output  file 
cond_out  (node_out) ; 
cond_out {sigma_out) ; 
cond_out {gamma_out) ; 

//write  vortex  filament  node  position  to  file 

out (x  n,y  n, z  n,node^out,max_node, n _ fil);  //main  vortex  tube 

out  ( r”x  n7t_y__n,  r_z_n ,  r_node_out ,  r__max_node ,  r__n_f  il ) ;  // radial  spokes 

//write  vortex  filament  core  radius  to  file 

out (sigma, sigma, sigma, sigma_out,max_node, n_fil) ;  / /main  vortex  tube 

out  (r_sigma,  r_sigma,  r_sigma,  r_sigma_out,  r_max_node,  r_n_fil) ;  //radial  spokes 

//write  vortex  filament  circulation  to  file 

out  (gamma,  gamma,  gamma,  gamma_out,max_node,  n_fil) ;  / /main  vortex  tube 

out  ( r_gamma ,  r_garoraa ,  r_gamma ,  r_gamma_out ,  r_max_node ,  r_n_f  i  1 ) ;  / /  radial  spokes 

f  ! - Begin  main  program  loop - 

while  (t_real  <  t_final-dt)  { 


calc  indO  ; 


//calculate  induced  velocity  on  vortex  filament  nodes 


t_index  =  t_index  +  1; 
t  real  =  t  real  +  dt; 


//write  current  simulation  time  to  screen 
int  c2  =  clockO; 

cout  «  ”t_index:  "«  t_index«”  t_real:  "«t_real  « 
double {c2~cl) /CLK_TCK  «  endl; 


run  time:  "  « 


//write  output  to  file  .  ■  ^  r 

if  (t  index  >=  t_out_start  &&  t_index  ==  out__count*out_sJcip)  { 

cout  «endl;  .  ^  ^  . 

cout«"Writing  node  position  to  file,  t_index=  "«  t_index«endl; 

cout  «endl; 


//write  vortex  filament  node  position  to  file 

out(x  n,y  n,z  n, node_out,raax_node, n_fil) ;  //main  vortex  tube 

out ( r“x_n, r_y_n, r_z_n, r__node_out , r_max_node , r_n_f il ) ;  //radial  spokes 


//write  vortex  filament  core  radius  to  file 

out (sigma, sigma, sigma, sigma_out,max_node, n_fil) ;  / /main  vortex  tube 
out (r_sigma, r_sigma, r_sigma, r__sigma_out, r_max_node, r_n_fil) ;  //radial 


spokes 


//write  vortex  filament  circulation  to  file 

out  (gamma,  gamma, gamma,  gamma_out,max_node,  n_fil) ;  //main  vortex  tube 

out(r  gamma,  r_gainma,r_gainma,r_gairana_out,r_max_node,r_n_fil);  //radial  spokes 


out_count  =  out_count  +  1;  //index  for  writing  output  to  file 

) 

) 

int  c3  =  clockO; 
cout  «  endl; 

cout  «"Total  run  time:  ”«  double (c3“Cl) /CLK_TCK  «  endl; 
cout  «  endl  «  "Press  any  key  to  continue..."; 
getch ( ) ; 


outf ile . close ( ) ; 


return  0; 

} 

// - 


/ 1 _ Function  read_data - 

//Read  data  from  vlinec.inp.  Input  file  requires  '$'  characters  and  end  of  line 

- - 

int  read_data  ( ) 

^  ifstream  infile ("vlinec_rad_spoke .inp",  ios::in); 

if  (iinfile)  return  0; 
infile .ignore (10000, 

infile  »  in_case; infile .ignore (1000,  *\n*); 

infile  »  t_final; infile . ignore (1000,  '\n'); 

infile  »  U_inf;  inf  ile.  ignore  (1000,  '\nM; 
infile  »  dz; inf ile. ignore (1000,  '\n'); 

infile  »  dt; infile .ignore (1000,  '\n'); 

infile  »  max_node; inf ile. ignore (1000,  ‘\n’); 

infile  »  r_fil_in; inf ile. ignore (1000,  *\n’); 

infile  »  n_rad; inf ile. ignore (1000,  ’\n'); 

infile  »  t_out_start;infile. ignore (1000,  '\nM; 

infile  »  out_skip; inf ile. ignore (1000,  *\n’); 

infile  »  r_tube0; inf ile. ignore (1000,  '\n'); 

infile  »  r_tubel; infile .ignore (1000,  '\n*); 

infile  »  r_tube2; inf ile. ignore (1000,  ’\n'); 

infile  »  r_tube3; inf ile. ignore (1000,  ’\n’); 


infile 

infile 

infile 

infile 

infile 

infile 

infile 

infile 

infile 

infile 

infile 

infile 

infile 

infile 

infile 

infile 

infile 

infile 

infile 

infile 

return 


»  r  tube4 ; inf ile . ignore ( 1000 ,  ^\n’)/ 

»  r  tubeS; inf ile . ignore {1000/  '\n')f 

»  r__tube6;  inf  ile.  ignore  (1000/  '\nM; 

»  r^tubeT; inf ile. ignore (1000,  ’\n'); 

»  r^ubeS; infile. ignore  (1000/  ’\nM/ 

»  r  tube9; inf ile, ignore (1000/  *\n')/ 

»  r~tubel0; inf ile. ignore (1000/  ’\n'); 

»  r_tubell; infile. ignore (1000/  *\n'); 

»  r_tubel2; inf ile. ignore (1000/  '\n'); 

»  sigina_init;  infile .  ignore  (1000/  '\n')/' 

»  C;infile. ignore  (1000,  *\nM; 

»  start_C_inc; infile . ignore ( 1000 , 

»  stop_C_inc; infile. ignore (1000/ 

»  c_inin;  inf  ile.  ignore  (1000,  '\n'); 

»  alpha;infile. ignore  (1000/  '\nM; 

»  z_on;  inf  ile.  ignore  (1000/  ’\nM/ 

»  z”init;infile. ignore  (1000,  '\nM; 

»  x_init;  inf  ile.  ignore  (1000,  '\nM; 

»  y^init; infile .ignore (1000,  '\n*); 

. close ( ) ; 

0; 


*\nM; 

\n’); 


} 


II - Function  init 

//Initialize  variables 

// - 


int  init  ( ) 

r_fil  =  new  int  [n_rad] ; 
r  tube  =  new  double  [n_rad] ; 
n^fil  =  0/ 


//initialize  number  of  filaments  in  a  given  radius 


r_fil[0]  - 
r_fil[l]  = 
r_fil[2]  = 
r_fil[3]  - 
r_fil[4]  = 
r_fil[5]  = 
r_fil[6]  - 
r_fil[7]  - 
r_fil[8]  = 
//r  fil[9] 


r_fil_in; 
r_fil_in; 
r_fil_in; 
r_fil_in; 
r_fil_in; 
r_fil_iii; 
r_fil_in; 
r_fil_in; 
r_fil_in; 
r  fil  in; 


//r_fil[10]=  r_fil_in; 

//r_fil[ll]=  r_fil_in; 

//r_fil[12]-  r__fil_in; 

//Initiallize  vortex  filament  radii 

r_tube [0] ^r_tube0; 

r_tube  [  1  ]  ==r_tube  1  ; 

r_tube [ 2 ] =r_tube2 ; 

r_tube [3] =r_tube3 ; 

r_tube  [  4  ]  =r_tube4  ; 

r_tube  [5]  =r_tube5 ; 

r_tube  [  6]  =r_tube  6 ; 

r^ube  [7]  =r_tube7 ; 

r_tube [ 8 ] =r_tube8 / 

/ / r _t ub  e [ 9 1 =  r_t ube  9 ; 

/ / r_tube [10] =r_t ube 1 0 ; 

/ /r_tube [11] =r_tube 1 1 ; 

/ / r_t ube [12] =r_tube 1 2 ; 

for  (int  j  -  0;  j  <=  n^rad-l;  j++) { 
n_fil  =  r_fil[j]+n_fil; 

} 


n_fil  -  n__fil  +  1;  //add  centerline  filament 
double  **seg_length  =  new  double  *[max_node]; 

for  (int  i  =  O;i<max_node;i++) {seg_length [i]^  new  double  [n_fil];} 

//initialize  simulation  time 
t_real  -  0.0; 
t_index  =0/ 
out_count  =  1; 

pi  =  M__PI; //define  pi 

2_init  ==  z_init  +  dz; 

//Initiallize  dynamic  arrays 
core__con  -  new  double  [n_fil]  ; 

on_node  ^  new  int  [n_fil] ; 

gamma  =  new  double  *[max_node]; 

for  {int  i  =  0;i<max_node; i++) {gamma [i] =  new  double  [n_fil];} 
x_n  =  new  double  *[max_node]; 

for  {int  i  =  0;i<max_node;i++)  {x_n  [i]  ==  new  double  [n__fil] ; } 
y_n  =  new  double  *  [max__node]  ; 

for  (int  i  =  0/i<max_node; i++) {y_n[i] =  new  double  [n_fil];} 
z_n  =  new  double  *[max_node]; 

for  {int  i  -  0;i<max__node;i++)  {2_n[i]=  new  double  [n_fil];} 

u_induced_total_x  =  new  double  * [max_node] ; 

for  (int  i  =  0;i<max_node; i++) {u_induced_total_x [i] = 

new  double  [n_fil];} 

u_induced_total_y  =  new  double  * [max_node] ; 

for  (int  i  =  0;i<max_node;i++) {u_induced_total_y [i]= 

new  double  [n^fil];} 

u_induced_total_z  =  new  double  *[max_node3; 

for  (int  i  =  0;i<max_node;i++) {u_induced_total_z [i] = 

new  double  [n_fil];} 

sigma  =  new  double  *[max_node]; 

for  (int  i  =  O;i<max_node;i++) {sigma [i] =  new  double  [n^fil];} 

//radial  spokes  arrays  and  locations 

int  m_2  =  10; 
n_spokes  =  16; 
r_n_fil  -  10*n_spokes; 
r_max_node  =  5/0.1; 

r_gamma  =  new  double  *  [  r_max_node  ]  ; 

for  (int  i  =  0; i<r_max_node; i++)  {r_gainma [i]  =  new  double  [r_n_fil];} 
r_x_n  =  new  double  * [r_max_node] ; 

for  (int  i  =  O;i<r_max_node;i++)  {r_x__n[i3=  new  double  [r_n_fil];} 
r_y_n  ^  new  double  * [r_max_node] ; 

for  (int  i  =  0; i<r_max_node;i++) {r_y_n [i] =  new  double  [r_n_fil];} 
r_z_n  =  new  double  * [r_max_node] ; 

for  (int  i  =  0; i<r_max_node; i++) {r_z_n [i] =  new  double  [r_n_fil];} 

r_u_induced_total_x  =  new  double  *  [r__max_node]  ; 

for  (int  i  =  0; i<r_max_node; i++)  {r__u_induced_total_x [i]  = 
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new  double  [r_n_fil];} 

r  u_induced_total_y  -  new  double  *  [rjnax_node]  ; 
for~{int  i  =  O;i<r_inax_node;i++)  {r_u_induced_total_y[i]  = 
new  double  [r_n__fil];} 

r  u  induced_total_z  =  new  double  *  [niax__node]  t 

for”'{int  i  =  O;i<rjnaax_node;i++)  {r_u_induced_total_z  [i]  = 

new  double  [r_n_fil];} 
r  sigma  -  new  double  *  [r_max_node] ; 

for  (int  i  =  O;i<r_max__node/ i++)  {3^_sigma [i]  =  new  double  [r_n_fxlj,} 
r  core  con  -  new  double  *  [r_max_node] ; 

for  (int  i  =  0; i<r_max_node; i++) {r_core_con[i] =  new  double  [r_n_fil] / } 

r  seg  length  =  new  double  *  [ r_iiiax_node ] ;  \ 

for  (int  i  =  0; i<r_jnax_node;i++) {r_seg_length [i] -  new  double  [r_n_fxl] ; } 

//Initialize  x_n,  y_n,  z_n 

//Note:  z_n  is  to  the  right  (axial  direction),  x_n  xs  up, 

//y_n  is  out  of  page.  When  plotting,  plot  z_n  to  the  right,  x_n  up, 
//and  y_n  out  of  page. 

//subroutine  in  init_incase_rad_spoke . cpp  that  initializes  vortex  fxlament 
//node  location 

init_incase_rad_spoke  (in_case,  x_n,  y_n,  z_n,  gamma,  C, max_node ,  r_f xl, n_f il,  n_raa, 
r  tube,x  init,  y_init,  z_init,dz,  start_C__inc,  stop_C_inc,  C_min, 

r  gamma ,  r_x_n ,  r_y_n ,  r__z__n ,  r_max_node ,  n_spokes ,  m_z ) ; 


//Intialize  sigma,  find  core  constant  to  maintain  constant  volume  in 
//vortex  filament  segment.  This  allows  simulation  to  handle  vortex 
//stretching 

for  (int  n  =  0;  n  <=  n_fil”l;  n++) { 

for  (int  m  =  0;  m  <=  max_node-2;  m++)  { 

sigma  [m]  [n]  =  sigma_init; 

seg_length[m] [n]  -  sqrt ( (x_n[m+l] [n] -x_n[m] [n] ) * 

(X  n[m+l]  [n]-x_n[m]  [n] )  +  (y_n  [m+1]  [n]-y_n[m]  [n} )  *  (y_n  [m+1]  [n] - 
y_n[m]  [n] )  +  (z_n[m+l]  [n]  '-z_n[m]  [n] )  *  (z_n[m+l]  [n]  -z_n[m]  [n] )  )  ; 

} 

core_con[n]  =  sigma [1] [n] *sigma[l] [n] *seg_length[l] [n] ; 

} 

//intialize  sigma  for  radial  spokes 
for  (int  n  =  0;  n  <=^  r__n_fil-l;  n++)  { 

for  (int  m  -  0;  m  <=  r_max_node-'2;  m++)  { 
r_sigma[m]  [n]  =  sigma_init; 

r  seg_^length [m]  [n]  =  sqrt(  (r_x_n[m+l]  [n]~r_x^n[m]  [n]  )* 

(7  X  n[m+l]  [n]-r_x_n[m]  [n] )  +  (r_y_n  [m+1]  [n]-r_y_n[m]  [n] )  * 
(r3y_n[m+l]  [n] -r__y_n [m]  [n] )  +  (r_z_n [m+1]  [n] -r_z__n[m]  [n] )  * 
(r_z_n[m+l)  [n] -r_z__n[m]  [n] ) )  ; 

r  core  con[m] [n]  =  r^sigma[m] [n] *r_sigma [m] [n] * 
r_seg_length [m] [n] ; 

} 

} 

return  0; 


} 
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j! - Function  cond_out -  “ 

//  Write  inputs  to  1st  row  of  output  file 

- - 

void  cond_out ( char  *out) 

of stream  out file (out,  ios::app); 

//need  to  fill  34  columns 

//output  simulation  variable  to  top  row  of  output  file 

outfile  «  in  case  «  "  "«  n_fil  «  "  "«max_node«  "  "«dt  « 

«  dz  «  "  "«  0  inf  «  "  "«  t_out_start  «  "  "«  out_s)cip  «  " 

«  r  tubel  «  "  "«  r_tube2  «  "  "«  r_tube3  «  "  ”  «r_tube4« 

«  r"tube5  «  "  "«  sigma  init  «  "  "«  C  «  "  "«  alpha  «  "  ” 

«  z'on  «  "  "«  n  rad  «  "  "«  r_fil[0]  «  "  "«  r_fil[13«  "  ” 

«  r“fil[2]«  "  '•  «r  filt3]«  "  "«r_fil[4]«  "  "«  r_fil[5]«  ”  " 

«  r"fil[6]«  "  "  «r_fil[7]«  "  "«r_fil[8]«  "  "«  r_fil[9]«  " 

for  (int  p=29  ;  p<=n_fil+l  ;  p++) {outfile  «  0  «  " 
outfile«  endl; 


} 

f  j - Function  out -  "* 

//  Write  x_n,  y_n,  z_n  to  file 

- - 

int  out(double**a,  double**b,  double**c,  char*out,  int  m_row,  int  n_col) 

( 

ofstream  outfile(out,  ios::app)/ 
if  (! outfile)  return  0; 

//output  x_n 

for  (int  m  ==  0;  m  <=  m_row-l;  m++)  { 

for  (int  n  ==  0;  n  <==  n_col-l;  n++)  { 
out file. precision (6) ; 
outfile  «  a[m]  [n]  «  " 

} 

outfile  «  t_index  «endl; 

} 

//output  y_n 

for  (int  m  =  0;  m  <=  m_row-l;  m++) { 

for  (int  n  -  0/  n  <-  n_col-l;  n++) { 
outfile. precision (6) ; 
outfile  «  b[m][n]  «  ” 

} 

outfile  «  t_index  «endl; 

} 

//output  z_n 

for  (int  m  =  0/  m  <=  m_row-l;  m++) { 

for  (int  n  =  0;  n  <=  n_col-l;  n++) { 
outfile. precision (6) ; 
outfile  «  c[m][n]  «  ” 

} 

outfile  «  t_index  «endl; 

} 

//outfile . close { ) ; 
return  0; 

} 

// - Function  calc_ind - 

//  Perform  induced  velocity  calculations 

U - - - 

void  calc_ind ( ) 

{ 

init_ind ( ) ; 

//calculate  induced  velocity  on  each  node  from  every  other  vortex 
//filament  segment 
int  n; 

for  (int  jl  =  0;  jl  <=  n_rad;  jl++){ 
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if  {jl==0)  {n  =  0;} 

else  if  (jl==l){n  =  1;} 

else  if  (jl  >  l){n  =  r_fil[jl-2]+  n;} 

for  (int  m  =  0;  m  <=  on_node [n] -1;  m++) { 
for  (int  k  =  0;  k  <=  n_fil-l;  k++) { 

for  (int  j  =  0;  j  <=  on_node [ k] -2 /  j++) { 


double 

double 

double 

double 

double 

double 

double 

double 

double 

double 

double 

double 


double 

double 

double 


gain_eq/ 

dlx  =  x_nEj+l]  [k]-x__n[j]  [k]; 
dly  =  y_n[j+l]  [k]-y__n[j]  [k]  ; 
dlz  =  z_n[j+l] [k]-z_n[j] [k] ; 

rx  =  X  n[m][n]  -  (0.5*  (x_n[j+l]  [k]  -  x_n[j][k3)  +  x_.n[j][k]); 

ry  =  yn[m][n]  -  (0.5*  (y_n[j+l]  [k]  -y__n[j][k3)  +y_n[j][k]); 

rz  =  z_n[m}[n]  -  (0 . 5*  (z_n[  j+1]  [k]  -  z_n[j][k])  +  z_n[j][k]); 

r_cross_dr_x  ==  rz*dly  -  ry*dlz; 

r_cross_dr_y  =  -(rz*dlx  -  rx*dlz) ; 
r_cross_dr_z  =  ry*dlx  -  rx*dly; 

abs_rk  =  sqrt(rx*rx  +  ry*ry  +  rz*rz) ; 

s  =  l/(pow(l+ (alpha*  (sigma  [m]  [n]  *sigina[m]  [n]  +  sigma  [j  ]  [k]  *sigma  [j  ]  [k] )  / 

(2*abs_rk*abs_rk) ) , 1.5) ) ; 

if  (0.5*  (z_n[j+l]  tk]  -  z_n[j][k3)  +  z_n[j]  [k]  <=  0)  { 
gam_eq  =  0.0; 

} 

else  { 

gam__eq  =  gamma  [  j  ]  [  k]  ; 

} 

u_ind_dummy_x  =  (-gam_eq/ (4 . 0*pi) ) *r_cross_dr_x*s/ (pow(abs_rk, 3) ) ; 
u_ind_dummy_y  =  (-gam_eq/ (4 . 0*pi) ) *r_cross_dr_y*s/ (pow(abs_rk, 3) ) ; 
u_ind_duinmy_z  =  {~gam_eq/ (4 . 0*pi)  )  *r_cross_dr_z*s/ (pow(abs_rk,  3)  )  ; 

if  (j  ==  0  &&  k  -=  0)  { 

u_induced__x  -  u_ind_dummy_x  ; 
u_induced_y  =  u_ind_duirany_y; 
u  induced__z  =  u_ind_duinmy_z ; 

else  { 

u  induced_x  =  u_induced__x  +  u_ind_dummy__x; 
u  induced_y  =  u_induced__y  +  u_ind_duinmy_y; 
u  induced  z  =  u_^induced_z  +  u_ind_dummy_z ; 

}“ 

} 

} 

if  (z_n[m]  [n]  <  0)  { 

u_induced_total_x [m] [n]  -  0; 
u_induced_total_y [m] [n]  =  0; 
u_induced__total_z  [m]  [n]  =  0; 

} 


else  { 

u_induced_total_x [m] [n] 
u_induced_total_y  [m]  [n] 
u_induced_total_z [m] [n] 
} 


u  induced_total_x [m] [n]  +  u_induced_x; 
u”induced”total_y[m3  [n]  +  u__induced_y; 
u  induced  total  z [m] [n]  +  u_induced_z; 


197 


) 

void  r 
{ 


double 

double 

double 

double 

double 

double 

double 

double 

double 

double 

double 

double 


} 

} 

//calculate  induced  velocities  on  and  from  radial  spoke  calculation 
r_calc_ind  ( ) ; 

//calculate  induced  velocity  for  all  nodes,  based  on  symmetry 
sym_calc { ) ; 
r_sym_calc ( ) ; 
u_update  ( ) ; 


calc_ind { ) 

u_induced_x  =0; 
u_induced_y  =  0; 
u_induced_z  =  0; 

r_on_fil  -  0; 

for  (int  n  =  0;  n  <=  r_n_fil‘-l;  n++)  { 
if  (r_2_n[0] [n]  >=  0)  { 

r_on_fil  =  r_on_fil  +  1  ; 

} 

} 

if  (r_on_fil  >=  r_n_fil) { 
r_on_fil  =  r_n_fil; 

} 

//calculate  velocity  induced  on  vortex  tube  filaments  by  spokes 
int  n; 

for  (int  jl  =  0;  jl  <=  n_rad;  jl++) { 

if  (jl==0) {n  =  0;} 

else  if  (jl-=l){n  =  1;} 

else  if  (jl  >  l){n  =  r_fil[jl-2]+  n; } 

for  (int  m  =  0;  m  <=  on__node [n] -1;  m++)  { 

for  (int  k  =  0;  k  <=  r_on_fil’-l;  k++)  { 

for  (int  j  =  0;  j  <==  r_max_node-2;  j++)  { 

gam_eq; 

dlx  =  r_x_n[j+l]  [k]“r_x_n[j]  [k]  ; 
dly  =  r_y_n[j+l] [k]-r_y_n[j] [k] ; 
dlz  =  r_z_n[j+l} [k]-r_z_n[j3  [k] ; 

rx  =  X  n[m][n]  -  (0.5*  (r_x_n[j+l]  [k]  -  r_x__n[ j]  [k] )  +  r_x_n[ j ]  [k] )  ; 

ry  =  y"n[m][n]  -  (0.5* {r_y_n[ j+1] [k]  -  r_y_n[j][k])  +  r_y_n [ j ] [k] ) ; 

rz  =  z_n[m][n]  -  (0 . 5* (r_z_n[ j+1] [k]  -  r_z_n[j3[k])  +  r_z_n [ j ] [k] ) ; 

r_cross_dr_x  =  rz*dly  -  ry*dlz; 
r_cross_dr_y  =  “(rz*dlx  -  rx*dlz) ; 
r_cross_dr_z  =  ry*dlx  -  rx*dly; 

abs_rk  =  sqrt(rx*rx  +  ry*ry  +  rz*rz); 

s  =  l/(pow(l+(alpha*(sigma[m] [n3*sigma[m] [n]  +  r_sigma [ j ] [k] * 

r_sigma[j3  [k] )  /  (2*abs__rk*abs_rk) )  ,1.5)  ) ; 

if  (0.5*  (r_z_n[ j+1]  [k]  -  r_2_n[j3Ck3)  +  r_z_n[j][k]  <=  0)  { 
gam_eq  =  0.0; 

} 

else  { 

gam_eq  =  “r_gainma[j]  [k]; 
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} 


double  u_ind_duinmy_x  == 
double  u_ind_duinmy_y  = 
double  u_ind_duinmy_z  = 


(-gain_eq/ (4 .0*pi) )  *r_cross_dr_x*s/ (pow(abs_rk;,  3) ) ; 
{-gam~eq/  (4 .0*pi) )  *r_cross_dr_y*s/  (pow (abs__rk,  3) ) ; 

(-gam_eq/(4.0*pi) )*r_cross_dr_z*s/(pow(abs_rk,3) ) ; 


if  (j  ==  0  &&  k  ==  0)  { 


u_induced_x  =  u_ind_duinmy_x; 
u_induced_y  =  u_ind_duitimy_y ; 
u  induced  z  =  u_^ind_duinniy_z  / 

}” 

else  { 

u  induced_x  =  u_induced_x  +  u_ind_dununy_x; 
u”~induced_y  -  u_induced_y  +  u_ind__duinmy_y/ 
u~induced_z  =  u_induced_z  +  u_ind_duinmy_z  ; 

}" 

} 

} 

if  (z_n[m] [n]  <  0) { 

u_induced_total__x [m]  [n]  =  0; 
u_induced_total_y [m] [n]  =  0; 
u_induced_total_z [ml [n]  =  0; 

} 

else  { 

u_induced_total__x  [m]  [n] 
u_induced_total_y  [m]  [n] 
u_induced_total_z [m] [n] 

} 

} 

} 

//calculate  velocity  induced  on  spokes  by  vortex  tube  filaments 


=  u_induced_total_x [m] [n] 
=  u_induced_total__y  [m]  [n] 
=  u_induced_total_z [m] [n] 


+  u__induced_x; 
+  u_induced_y; 
+  u  induced  z; 


for  (int  n  =  0;  n  <=  r_n_fil“l;  n++) { 

for  (int  m  =  0;  m  <=  r_max_node-l;  m++)  { 
r_u_induced_total_x  [m]  [n]  =  0; 
r_u_induced_total_y [m] [n]  -  0; 
r_u_induced_total_z [ml [nl  =  0; 

} 

} 

u_induced__x  =  0; 
u_induced_y  =0; 
u  induced_z  =0; 


for  (int  jl  =  0;  jl  <=  r_on_fil/n_spokes-l;  jl++){ 
n  =  n_spokes*jl; 

for  (int  m  =  0;  m  <=  r_max_node“l;  m++)  { 

for  (int  k  =  0;  k  <=  n_fil-l;  k++) { 

for  (int  j  =  0;  j  <=  on_node [k] -2;  j++) { 


double  gam__eq; 

double  dlx  =  x__n  [  j+l]  [k] -x_n  [  j  ]  [k]  ; 
double  dly  =  y_n [ j+1] [k] -y_n[ j ] [k] ; 
double  dlz  =  z_n[ j+1] [k] -z_n[ j] [k] ; 

=  r  X  n[m]  [n]  -  (0.5*  (x_n[ j+11  [k]  -  x_n[j][kl)  +  x_n[j][k]). 


double  rx 


double  ry  =  r_y_ntm]  [n]  -  (0.5*  (y_n[j+l]  [k]  -  y_n[j][k])  +  y_n[j][k]); 
double  rz  =  r_z_n[m]  [n]  -  (0.5*  (z_n[j+l]  [k]  -  z_n(j][k])  +  z_n[j][k]); 

double  r_cross__dr__x  =  rz*dly  -  ry*dlz; 
double  r__cross_dr_y  =  -(rz^dlx  -  rx*dlz); 
double  r_cross_dr_z  =  ry*dlx  -  rx*dly; 

double  abs_rk  =  sqrt(rx*rx  +  ry*ry  +  rz*rz); 

double  s  =  1/ (pow(l+ (alpha* {r_sigma[m] [n]*r_sigma[m] [n]  + 

sigma [j] [k]*sigma[j] [k] ) / {2*abs_rk*abs_rk) ) ,1,5) ); 

if  (0.5*{z_n[j+l]  [k]  -*  z_n[j][k])  +  z_n[j]  [k]  <=  0)  { 
gain_eq  =  0.0; 

} 

else  { 

gam_eq  ==  gainma[j]  [k]  ; 

} 

double  u  ind  duinmy_x  =  (-gam_eq/  (4.0*pi) )  *r__cross_dr_x*s/ {pow(abs_rk,  3)  )  , 
double  u_ind_duinmy_y  =  (-gam_eq/  (4 .0*pi) )  *r_cross_dr__y*s/  {pow{abs_rk,  3) ) ; 
double  u_ind__dummy_z  =  (-~gam_eq/  (4 .0*pi) )  *r_cross_dr_z*s/  (pow{abs_rk,  3) ) ; 

if  (j  ==  0  &&  k  ==-  0)  { 

u_induced_x  =  u_ind_duinmy_x; 
u_induced_y  =  u_ind_duinmy__y ; 
u  induced__z  =  u__ind_dummy_z ; 

}” 

else  { 

u_induced_x  =  u_induced_x  +  u_ind_duinmy_x ; 
u  induced__y  =  u_induced_y  +  u_ind_duinmy__y ; 
u  induced  z  =  u_induced__z  +  u_ind_duinmy_z ; 

}” 

} 

} 

if  (r_z_n[in]  [n]  <  0)  { 

r_u_induced_total_x [m] [n]  =0; 
r_u_induced_total_y  [m]  [n]  =  0; 
r_u_induced_total_z [m] [n]  =0; 

} 

else  { 

r_u_induced_total_x  [m]  [n] 
r_u_induced__total_y  [m]  [n] 
r_u_induced_total_z [m] [n] 

} 

} 

} 

} 

// - Function  init_ind - 

//  Initialize  induced  velocity  calculations 

- - 

void  init_ind() 

{ 

for  (int  n  =  0;  n  <=  n_fil“l;  n++) { 
on_node[n3  =  0; 

for  (int  m  =  0;  m  <-  xnax_node-l;  m++)  { 


u_induced_total_x [m] [n]  =  0; 
u  induce d_total_y [m] [n]  -  0; 


=  r  u_induced_total_x  [m]  [n] +u__induced__x; 
~  r__u_induced_total_y  [m]  [n] +u_induced_y; 
=  r  u  induced  total_z[m]  [n] +u_induced__z/ 


u_induced_total_z  [iti]  [n]  -  0; 

if  {z_n[m] [n]  >=  0)  { 

on_node[n]  =  on_node[n]  +  1; 

} 


} 


if 

} 


(on_node[n]  >-  max_node>  { 
on  node[n]  -  max_node; 


11 - Function  init_ind -  " 

//  Update  node  positions  and  filament  parameters 

- - 

void  u_update ( ) { 

double  **seg_length  =  new  double  *[max_node]; 

for  (int  i  =  0;i<max_node;i++) {seg_lengthti] =  new  double  [n_fil];} 


for  (int  n  =  0;  n  <=  n_fil“l;  n++) { 

for  (int  m  =  0;  m  <=  max_node-l;  m++)  { 

X  n[m]  [n]  ==  x_n[m]  [n]  +  u_induced_total_x  [m]  [n]  *dt; 
y  n[m][n]  ^  y^n[m] [n]  +  u_induced_total_y [m] [n] *dt; 
z  n[m] [n]  =  z_n[m3 [n]  +  (u_induced_total_z [m] [n]  +  U_inf)*dt; 

}” 

} 

for  (int  n  =  0;  n  <=  n_fil“l;  n++) { 

for  (int  m  =  0;  m  <=  max_node“2;  m++) { 

seg  length[m][n]  =  sqrt  ( (x_n[m+13  [n]--x_n[m]  [n] )  * 

(x_n[m+l]  [n]-x_n[m]  [n] )  +  (y__n[m+13  [n]-y_n[m3  [n] )  *  (y_n  [m+l]  [n] - 
y_n[m]  [n]  )  +  (z_n[m+13  Cn]-2_n[m]  [n]  )*(z_n[m+l]  [n]-z_n[m3  [n] ) )  ; 

sigma [m3 [n3  =  sqrt ( {core_con[n3 /seg_length[m3  [n3 ) ) ; 

} 

} 


//update  radial  spokes 
r  seg  length  =  new  double  *  [  r__max_node  3  ; 

for  (int  i  =  0;i<r_max__node;i++)  {r_seg_length[i3=  new  double  [r_n_fil3 ; } 

for  (int  n  =  0;  n  <=  r_n_fil-'l;  n++)  { 

for  (int  m  =  0;  m  <-  r__max_node-l;  m++)  { 

r_x_n[m3[n]  r_x_n[m3  [n3  +  r_u_induced_total_x  [m]  [n3  *dt; 
r  y  n[m3 [n]  =  r_y_n[m3  [n3  +  r_u_induced_total_y [m3 [n3*dt; 
r^z”n[m3[n3  ^  r_z_n[m3  [n3  +  (r_u_induced_total_z [m]  [n3  +  U_inf)*dt; 

} 

} 

for  (int  n  0;  n  <=  r_n_fil-l;  n++)  { 

for  (int  m  =  0;  m  <=  r_max_node-2;  m++) { 

r_seg_length  [m3  [n]  -  sqrt  ( (r_x_n  [m+13  [n3-r_x_n[m3  [n]  )* 
{r_x_n[m+13  [n3-r_x_n[m]  [n3 )  +  (r_y__n [m+13  [n3“r_y_n[m3  [n3 )  * 
(r_y_n[m+13  [n3-r_y_n[m]  (n3 )  +  (r_z_n [m+13  [n3  ■-'r_z_n[m3  [n3 )  * 
(r_z_n[m+13  [n3  -r_z_n[m3  [n3 ) ) ; 

r_sigma[m3  [n3  =  sqrt  ( (r_core_con[m3  [n3 /r_seg_length[m3  [n3 ) ) ; 

} 

} 


//cout  «"Umax:  "«Umax«"  Umin;  "«Umin«endl; 
delete []  seg_length; 
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void  sym_calc  { ) { 
int  n; 

for  (int  jl  -  1;  jl  <=  n_rad;  jl++){ 
if  (jl==l) {n  -  1;} 

else  if  (jl  >  l){n  r_fil(jl-2]+  n;} 

for  (int  k  =  n;  k  <=  n+r_fil [ j 1-1] -1;  k++) { 
for  (int  m  -  0;  m  <=  on_node [k] -1;  m++) { 

double  u_rad  =  (x_n[in]  [n] /sqrt  (x_n  [m]  [n]  *x_n[m]  Cn]+y_n[in3  [n]  * 

y  n  [m]  [n] ) )  *u_induced__total_x  [m]  [n]  +  (y__n  [m]  [n]  / sqrt  (x_n  [m]  [n]  * 
x""n[m]  [n3+y  n[ral  [n]*y_n[m]  [n] ) )  *u_induced_total_y  [m]  [n]  ; 
double  u^theta  =  -(y_n[m] [n] /sqrt (x_n [m] [n]*x_n[m] [n] +y_n[m] [n] * 
y  n[m] [n] > ) *u_induced_total_x [m] [n3+(x_n[m] [n] /sqrt (x_n [m] [n] * 
x“n[m3  [n]+y_n[ia]  tn]*y_n[in3  [n] )  )  *u_induced_total_y  [m]  [n] ; 
double  u_ax  =  u_induced_total_z [m] [n] ; 

u  induced  total_x[in3  [k]—  (x_n[ni]  [k] /sqrt  (x__n[iti]  [k]*x_n[m3  [k]  + 
y  n[in]  [k]  *y__n  [m]  [k] ) )  *  (u_^rad)  -  (y_n[iu]  [k] /sqrt  (x__n  [ift]  [k] 
x'nLm]  [k]+yln[m]  [k]*y_n[m]  [k] ) )  *  {u_theta)  ; 

u  induced_total_y  [m]  [k]=  (y_n[m3  [k] /sqrt  (x__n[m]  [k]*x_n[m]  [k]  + 

y“n[m]  [k]*y_n[m]  [k] )  )  *  (u_rad)  +  (x_n[m]  [k] /sqrt  (x_n[m]  [k]* 
x_nCm]  [k]+y_n[m]  Ck]*y_n[m]  [k]))*(u_theta)/ 
u_induced_total_z  [m]  [k]  =*=  n_ax; 

} 

} 

} 

} 

void  r_sym_calc 0 { 
int  n; 

for  (int  jl  ==  0;  jl  <=*  r_on_fil/n_spokes-l;  jl++){ 
n  =  n_spokes*jl; 

for  (int  k  =  n;  k  <=  n+n_spokes-l;  k++) { 
for  (int  m  =  0;  m  <=  r_max_node-l;  m+t) { 

double  u__rad  =  (r_x_n  [m]  [n] /sqrt  (r_x_n[m3  [n]  *r_x_n  [m]  [n] +r_y_n  [m]  [n]  * 

r  y  n  [m]  [n] ) )  *r_u__induced_total_x  [m]  [n]  +  (r_y_n  [m]  [n]  / 
sqrt  (r_x_n  [m]  [n]  *r_x_n  [in]  [n]  +r_y_n  [m]  [n]  *r_y_n  [m]  [n] ) )  * 
r  u_induced_total_y  [m]  [n] ; 

double  u_theta  =  -(r_y_n[m]  [n]/sqrt(r_x_n[in]  [n]*r_x_n[m]  [n]  + 
r  y  n[m] tn]*r  y  n[m] [n] ) ) *r_u_induced_total_x[m] [n] + {r_x_n[m] [n] / 
sqrt(r_x_n[mlTnT*r_x_ntra]  [n]+r_y_n[in]  [n]*r_y_n[m]  [n]  ) )  * 
r_u_induced_total_y [m] [n] ; 

double  u_ax  =  r_u_induced_total_z [m] [n] ; 

r_u_induced_total_x[m]  [k]=  (r_x_n[m]  [k]/sqrt(r_x_n[m]  [k]* 

r  X  n[m]  [k]+r_y_n[m]  [k]*r_y_n[in]  [k]))*(u_rad)-(r_y_n[m]  [k]/ 

sqrt(r_x_n[m]  [k]  *r_x_n[in]  [k]+r_y_n[m]  tk]*r_y_n[m]  [k] ) )  *  (u_theta)  ; 
r  u_induced_total_y [m] [k]=  (r_y_n[m] [k] /sqrt (r_x_n[m] [k] * 

r  x  n[m]  (k]+r_y_n[m]  [k]  *r_y_n[m]  [k] ) )  *  (u_rad)  +  (r_x_n[m]  [k]  / 

sqrt(r_x_n[m]  [k]*r_x_ntm]  [k]+r_y_n[m]  [k]  *r_y_n[m]  [k] ) )  *  (u_theta) ; 
r_u_induced_total_z [m] [k]=  u_ax; 

} 

} 

} 

} 


8.2  initjncase.cpp  -  Vortex  filament  initialization  code 

#include  <vcl.h> 

#include  <f stream. h> 

#include  <istream.h> 

#include  <iostream.h> 

#iiiclude  <conio.h> 

#include  <fastmath.h> 

#pragma  hdrstop 

void  init_incase_rad_spoke(int  in_case, double**  x_n,  double**  y_n, 
double**  z_n, double**  gamma, double  C,int  max_node,int*  r_fil,int  n_£il, 
int  n  rad, double*  r_tube, double  x_init,  double  y_init, 

double  2_init, double  dz, double  start_C_inc, double  stop_C__inc, double  C_min, 
double**  r_gamma, double**  double**  double**  3:__z_n, 

int  r_max_node,  int  n_spokes,  int  m_z) 

{ 

double  pi  =  M_PI; 
double  Cl; 
double  C_up; 
double  C_dovm; 
double  z_up; 
double  z_down ; 
double  z_upl; 
double  z_downl; 

double  ampl; 
double  per; 
double  phase; 
double  slope; 
double  inter; 
double  gamsum; 
double  r_tubel; 
int  incr  =0; 

//r_tube__suin  =  r_tube_sm  -  r_tube[0]; 


//Switch  cases  for  C  distribution 
switch  (in_case)  { 

case  1 : 

//Initialize  vortex  filament  node  positions 
for  (int  m  =  0;  m  <=  max_node-l;  m++) { 
x_n[m] [0]=x_init; 
y__n[m]  [0]=y_init; 
z__n[m]  [0]=z__init  -  dz*  (m)  ; 
for  (int  j  =  0;  j  <=  n_rad”l;  j++) { 

for  (int  n  =  0;  n  <=  r_fil[j]-l;  n++) { 

X  n[ml [n+l+incr]  =  r_tube [ j ] *cos ( (n) *pi*2/r_fil [ j ] )  +  x_init 
y__nCm]  [n+l+incr]  =  r_tube  [  j  ]  *sin  ( (n)  *pi*2/r_fil  [  j  ]  )  +  y_init 
z_n[m] [n+l+incr]  =  z_init-dz* (m) ; 

} 

incr  =  incr  +  r_fil[j]; 


} 

incr  =0; 

//Initialize  circulation  in  vortex  tube,  must  be  done  in  conjunction 

//with  radial  spoke  circulation 

z_up  =  stop_C_inc; 

z__down  =  start_C_inc; 

z__upl  =  z_init-l; 

z_downl  =  z_init; 

C_up  =  C; 
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C_down  =  C_inin; 
incr  -  0; 

for  (int  m  =  0;  m  <=  max_node-l;  in++)  { 
if  (z_n[m]  [0]>=z_upl  &&  z_n  [m]  [0]  <-2_downl>  { 

Cl  =  {{C_down~0)/(z_downl-z_upl))*  {z_downl-z_n [m]  [0]  )  +  0; 

else  if  (z_n[m] [0]>=z_down  &&  z_n [m] [0] <z_upl) { 

Cl  =  C_down; 

llse  if  (z_n[in]  [0]>=z_up  &&  z_n [m]  [0] <z_down  ){ 

Cl  =  ( (C  up-C  down)  /  (z_down-z_up) )  *  (z__down-z_n[m]  [0] )  +  C_down; 

} 

else  if  {z_n[m] [0] <z_up) { 

Cl  =  C_up; 

} 

else{ 

Cl  =  0; 

} 

cout  «C1«  ” 

//Gaussian  vorticity  distribution,  calculated  explicitily  based  on  the 
//method  of  the  Lagrange  multiplier 
gamma [m]  [0]  =0; 

int  j  "  0;  .  .  r 

for  (int  i  =  0;  i  <=  r_fil[j]-l;  i++) { 
gamma  [m]  [i+ 1+0* r_fil  [j ]  ]  =  . 008114549976*C1/  //r_^0 

gamma [m]  [i+l+l*r_fil [j ] ]  =  . 004780605093*C1; 

gamma [m]  [i+l+2*r_^fil[j]}  =  . 01061059381*C1;  //r_2 

gamma [m] [i+l+3*r_fil[j]]  -  . 01097097189*Cl;  //r_3 

gamma [m] [i+ 1+4* r_fil [j ] ]  =  .01223904142*C1;  //r_4 

gamma [m] [i+l+5*r_fil[j]]  =  .009090954102*C1;  //r_5 

gamma [m] [i+l+6*r_fil [j ] ]  =  . 002429657471*C1;  //r_6 

gamma [m] [i+l+7*r_fil [j ] ]  =  . 002628594294*C1;  //r_7 

gamma [m] [i+l+8*r_fil[j] ]  =  . 001640452531*C1;  //r_8 

} 

} 

cout  «endl; 

//calculate  total  circulation 
gamsum  =  0; 

for  (int  n  =  0/  n  <=  n_fil;  n++) { 
gamsum  =  gamma  [max_node-l]  [n]  +  gamsiom; 

} 

//initialize  radial  spokes 

//radial  spoke  node  location  interior  to  vortex  tube 
incr  =0; 

for  (int  j  =  0;  j  <=  m_z  -1;  j++) { 

for  (int  n  =  0;  n  n_spokes-l;  n++) { 

for  (int  m  =  0;  m  <=  n_rad-l;  m++) { 

r_x_n[m] [n+incr]  =  r^tube [n_rad-l-m] *cos ( (n) *pi*2/n_spokes) 

+  x__init; 

r_y_n[m] [n+incr]  =  r_tube [n_rad-l-m3 *sin ( (n) *pi*2/n_spokes) 

+  y_init; 

r_z_n[m] [n+incr]  =  start_C_inC“dz* ( j ) ; 

} 

} 

incr  =  incr  +  n_spokes; 

//radial  spoke  node  location  exterior  to  vortex  tube 
incr  =  0; 

for  (int  j  =  0;  j  <=  m_2  -1;  j++) { 

for  (int  n  =  0;  n  <=  n_spokes-l;  n++) { 


for  (int  m  =  9;  m  <=  r_max_nocie-l;  m++)  { 

r  X  n[ra]  [n+incr]  =  0 . X* {m-4) *cos ( (n) *pi*2/n_spok€s)  +  x_xnit; 
r3y“n[m]  [n+incr]  -  0 . 1*  (m-4)  *sin  (  (n)  *pi*2/n_spokes)  + 
r_z_n[m]  [n+incr]  =  start_C_inc-ciz*  { j  )  ; 

} 


} 

incr 


incr  +  n_spokes; 


} 


//radial  spoke  circuaition  interior  to  vortex  tube 
incr  =  Of¬ 


fer 


{int 

for 


j  =  0;  j  <=  m_z  -1;  j++) { 

(int  n  =  0;  n  <=  n__spokes-l;  n++)  { 
r  gamma [0] [n+incr] =  {C-C_min) /lO* . 001640452531; 
r'gamma[l] [n+incr] =  (C-C_min) /lO* . 002 6285942 94+r_gamma [0] 
-  {C-C_min) /10*.002429657471+r_gamma[l] 
(C-C_min) /lO* . 009090954 102+r_gamma [2] 
(C”C_min) /lO* . 01223904 1420+r_gamma [3] 
(C-C__min)  /lO*  .010970971890+r_gamma[4] 
(C-C_min) /lO* . 010610593810+r_gamma [5] 
(C-C__min)  /lO*  .  004780605093+r_gainma  [6] 
(C“C_min)/10* .008114549976+r_gamma[7] 


r_gamma[2] [n+incr] - 
r_gamma[3] [n+incr] = 
r_gamma[4] [n+incr] = 
r_gamma[5] [n+incr] = 
r_gamma[6] [n+incr]  = 
r_gamma[7] [n+incr] = 
r_gamma[8]  [n+incr]  = 


[n+incr] ; 
[n+incr] ; 
[n+incr] ; 
[n+incr] ; 
[n+incr] ; 
[n+incr] ; 
[n+incr] ; 
[n+incr]  ; 


incr  =  incr  +  n_spokes; 


//radial  spoke  circuaition  exterior  to  vortex  tube 
incr  =  Of¬ 
fer  (int  j  =  0;  j  <=  m_z  -1;  j++) { 

for  (int  n  =  0;  n  n_spokes-l;  n++) { 

for  (int  m  =  9;  m  <=  r_max_node-l7  m++) { 

r_gamma[m] [n+incr] =  r_gamma[8] [n+incr] ; 

} 

} 

incr  =  incr  +  n_spokes; 

} 

break; 


case  2: 

//Initialize  vortex  filament  node  positions 
for  (int  m  =  0;  m  <=  max_node~l;  m++) { 
x_n[m]  [0]==x_init; 
y_n[m] [0]=y_init; 
z_n[m]  [0]=2_init  -  dz*  (m)  ; 
for  (int  j  =  0;  j  <=  n_rad“l;  j++) { 

for  (int  n  =  0;  n  <=  r_fil[j]-l;  n++) {  ^ 

X  n[m]  [n+l+incr]  =  r_tube  [  j]  *cos  (  (n)  *pi*2/r__fil[  j] )  +  x_:Lnxt 
y”n[m] [n+l+incr]  =  r_tube [ j ] *sin( (n) *pi*2/r_fil [ j ] )  +  y_init 
z_n[m]  [n+l+incr]  =  z_init-dz*  (m)  ; 

} 

incr  =  incr  +  r_fil[j]; 


} 

incr  =0; 

//Initialize  gamma 
z  up  =  stop_C_inc; 
z_down  =  start_C_inc; 
z_upl  =  z^init-l; 
z_downl  =  z_init; 

C_up  =  C; 

C__down  =  C_min; 

incr  =  Of¬ 
fer  (int  m  =  0;  m  <=  max_node-l;  m++)  { 

if  (z_n[m] [0]>=z_upl  &&  z_n[m] [0] <=z_downl) { 


Cl  =  ({c_down-0)/(z_downl-2_upl))*  {z_downl-z_n[m] [0] )  +  0; 

else  if  (z_n[in]  [0]>=z_down  &&  z_n[in]  [0] <z_upl)  { 

Cl  =  C_down; 

else  if  (z_n[m]  [0]>=z_up  &&  2_n[m]  [0]<z_down  ){ 

Cl  =  ((C_up-C_down)/{z_down-z_up))*  (z_down-z_n[m]  [0] )  +  C_down; 

} 

else  if  (z_n[m] [0]<z_up) { 

Cl  -  C_up; 

} 

else{ 

Cl  =  0; 

} 

cout  «C1«  " 


//Gaussian  vorticity  distribution,  calculated  explicitily. 
gamma  [m]  [  0  ]  =  0 ; 


int  :i  =  0; 

for  (int 

gamma [m]  [i+l+0*r__fil  [  j]  ] 
gamma  [m]  [i+l+l*r_fil  [  j  ]  ] 
gamma  [m]  [i+l+2*r_f  il  [  j  3  ] 
gamma  [m]  [i+l+3*r__fil  [  j  ]  ] 
gamma [m]  [i+l+4*r_fil  [ j  ]  ] 
gamma  [m]  [i+l+5*r_fil  [  j  3  3 
gamma  [m3  [i+l+6*r_fil  [  j  3  3 
gamma  [m3  [i+l+7*r_fil  [  j  3  3 
gamma  [m3  [i+l+8*r_fil  [  j  3  3 


i  =  0;  i  <=  r_fil[j3-l;  i++) { 
=  .008114549976*C1;  //r_0 

-  .004780605093*C1;  //r_l 

=  *01061059381*C1;  //r_2 

=  .01097097189*C1/  //r_3 

=  .01223904142*C1;  //r_4 

=  .009090954102*C1;  //r_5 

=  .002429657471*C1;  //r_6 

=  .002628594294*C1;  //r_7 

=  .001640452531*C1;  //r_8 


cout  «endl; 

//cout  «endl«"r_tube[03  =  ”«r_tube[03  «endl; 
//cout  «endl«”r_fil[03  =  "«r_fil[03  «endl; 
gamsum  =  0; 

for  (int  n  =  0;  n  <=  n_fil;  n++) { 
gamsum  -  gamma  [max_node- 13  [n]  +  gamsum; 

//cout  «endl«"n_fil  =  ”«n_fil  «endl; 
cout  «endl«”sum  of  gamma  =  "«gamsum  «endl; 


//initialize  radial  spokes 

//radial  spoke  node  location  interior  to  vortex  tube 
incr  =  0; 

for  (int  j  -  0;  j  <=  m_z  -1;  j++) { 

for  (int  n  =  0;  n  n_spokeS“l;  n++) { 

for  (int  m  =  0;  m  <-  n_rad-l;  m++) {  .  .  ,  x 

r  X  n[m3 [n+incr3  =  r_tube [n_rad~l-m3 *cos ( (n) *pi*2/n_spokes ) 

+  X  init; 

r  y  n[m3  [n+incr3  ==  r_tube  [n_rad-l-m3  *sin  (  (n)  *pi*2/n_spokes) 

+  y_init; 

r_z_n[m3 [n+incr3  =  start_C_inc-dz* ( j ) ; 

} 

} 

incr  =  incr  +  n_spokes; 

//radial  spoke  node  location  exterior  to  vortex  tube 
incr  =0; 

for  (int  j  =0;  j  <-  m_z  -1;  j++) { 

for  (int  n  0;  n  <==  n_spokes“l;  n++)  { 

for  (int  m  =  9/  m  <=  r_max_node-l;  m++) { 

r  X  n[m3[n+incr3  =  0 . 1* (m-4) *cos ( (n) *pi*2/n_spokes)  +  x_init 
r_y_n[m3 [n+incr3  =  0 . 1* (m-4 ) *sin ( (n) *pi*2/n_spokes)  +  y_init 
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r_z_n[m] [n+incr]  -  start_C_inc-dz* { j ) ; 

} 

} 

incr  =  incr  +  n_spokes; 

//radial  spoke  circualtion  interior  to  vortex  tube 
incr  =  0; 

for  (int  j  =  0;  j  <=  in_z  -1;  j++)  { 

for  (int  n  =  0;  n  <=  n_spokes-l;  n++) ( 

r  ganiiiia [0]  [n+incr]  =  (C-C_niin) /lO*. 001640452531; 

r“gamma[l]  [n+incr]=  (C-C_inin)  /lO* . 002628594294+r_gainma[0]  [n+incr]  ; 
r“gairana[2] [n+incr]=  (C-Cjmin) /10*.002429657471+r_gamma[l] [n+incr] ; 
r'"gairana[3]  [n+incr]  =  {C-C_inin)  /lO*  .009090954102+r_gairana  [2]  [n+incr]  ; 
r“gainma[4]  [n+incr]=  (C-C_inin)  /10*.012239041420+r_gartima[3]  [n+incr]  ; 
r'"gamraa[5]  [n+incr]=  (C-C_inin)  /10*.010970971890+r_gainma[4]  [n+incr] ; 
r'gainma[6]  [n+incr]=^  (C-Cjmin)  /10*.010610593810+r_ganima[5]  [n+incr]  ; 
r“ganmia[7]  [n+incr]=  (C-C_min)  /10*.004780605093+r_gairatia[6]  [n+incr]  ; 
rlgairana[8]  [n+incr]=  (C-C_min) /10*.008114549976+r_gaitma[7]  [n+incr]; 

} 

incr  =  incr  +  n_spokes; 

} 

//radial  spoke  circualtion  exterior  to  vortex  tube 
incr  =  0; 

for  {int  j  =  0;  3  <=  m_z  -1;  j++} { 

for  (int  n  =  0/  n  <=  n_spokes-l;  n++) { 

for  (int  m  =  9;  m  <=  r_niax_node“l;  in++)  { 

r_gairana [m]  [n+incr]  ==  r_gainma[8]  [n+incr]; 

} 

} 

incr  =  incr  +  n_spokes; 

} 

break; 

} 

} 
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8.3  mitjncase.cpp  -  Vortex  filament  input  file 


II  This  is  the  input  file  for  vlinec.exe  Do  not  modify  the  structure. 

//  To  change  an  input  value,  just  replace  the  existing  value 
//Do  not  delete  the  dollar  sign  symbols  below 
- - 

$  ... 

1  //Simulation  initial  condition  case,  see  inxt_rncase . cpp 

2.7  //t_final.  .  .  final  simulation  time 

5  //U  inf.  .  .freestream  velocity 

0.1  //dz,  .  .  initial  axial  separation  between  vortex  filament  nodes 

0.01  //dt.  .  .simulation  time  step 

251  //max_node.  .  maxnxamber  of  nodes  in  a  filament  (excluding  radial  spokes) 
16  //numer  of  filaments  at  a  given  radius 
9  //number  of  radii  (excluding  centerline) 

0  //t  out  start.  .  .  t^index  output  starts 

30  //out  skip.  .  .  number  of  iterations  before  writing  data  to  file 

0.42  //r  tubeO.  .  .  initial  radial  location  of  an  axial  vortex  filament 

0.36  //r_tubel 

0.3  //r_tube2 

0.24  //r_tube3 

0.18  //r_tube4 

0.12  //r_tube5 

0.08  //r_tube6 

0.06  //r_tube7 

0.03  //r_tube8 

0.0  //r_tube9 

0.0  //r_tubel0 

0.0  //r_tubell 

0.0  //r_tubel2 

0.12  //sigma  init.  .  .  initial  vortex  filament  radius 
15.0  //Cmax.  .  .  final  value  of  total  circulation  in  vortex  tube 
-3  //start  C  inc.  .  .  axial  location  where  linear  C  increase  starts 

-4  //stop  C  inc  .  .  .  axial  location  where  linear  C  increase  ends 

0.1  //C  min.  .  .  intial  value  of  total  circulation  in  vortex  tube 

0.413  //alpha.  .  .  parameter  defining  velocity  near  filament  center 
0.0  //z  on,  .  .  axial  location  where  self  induction  starts 

0.0  //z"~init.  .  .  sets  overvall  initial  location  of  vortex  tube 

0.0  //x'^init.  .  .  sets  overvall  initial  location  of  vortex  tube 

0.0  //y  init .  .  .  sets  overvall  initial  location  of  vortex  tube 


9.  Appendix  B— Additional  Vortex  filament  simulation  results 
9.1  Vortex  filament,  r  =  0.12,  side  view 


Vortex  Filaments,  Time:  0.3 


2  3  4  5  6 

Vortex  Filaments,  Time:  1 .5 


Vortex  Filaments,  Time:  2.1 


Filaments,  Time:  2.7 
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9.2  Vortex  filament,  r  =  0.12,  rendered  perspective  view 

Vortex  Filaments,  Time:  0 


10 


-1 
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Vortex  Filaments,  Time:  1 .2 


10 


-1 


Vortex  Filaments,  Time:  1.8 


Vortex  Filaments,  Time:  2.1 
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Vortex  Filaments,  Time:  2.4 


10 


-1 


9.3  Vortex  filament,  r  =  0.12,  closer  rendered  perspective  view 

Vottt«  narnsns.  Tima;  0.  tWt  0.12  Vortax  FUamanta,  Tima:  O.S.  iWt:  0.12 


Vortex  Filamenfs,  Ttme;  2.4,  rir^:  0.12 


Vortex  RIameflte,  Time:  2.7,  rinit  0.12 


9.5  Color  contour  axial  velocity  plots 

Axial  Velocity,  meridonal  plane,  Time:  0.3 


Axial  Velocity,  meridonal  plane.  Time:  0.9 
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9.6  Color  contour  radial  velocity  plots 


Radial  Velocity,  meiidonal  plane,  Time:  0.3 
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Radial  Velocity,  meridonal  plane,  Time:  2.1 
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Azimuthal  Vorticity,  meridonal  plane,  Time:  2.1 


9.9  Color  contour  radial  vorticity  plots 

Radial  Vorticity,  meridonal  plane,  Time;  0.3 
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Radial  Vorticity,  meridonal  plane.  Time:  1.5 


